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ABSTRACT 

The  Electrical  and  Computer  Engineering  Department  at  the  Naval  Postgrad- 
uate School  has  a  need  for  additional  software  to  be  used  in  instructing  students 
studying  signal  processing.  This  software  will  be  used  in  a  PC  lab  or  at  home. 

This  thesis  provides  a  set  of  disks  written  in  APL  (A  Programming  Language) 
which  allows  the  user  to  input  arbitrary  signals  from  a  disk,  to  perform  various 
signal  processing  operations,  to  plot  the  results,  and  to  save  them  without  the  need 
for  complicated  programming. 

The  software  is  in  the  form  of  a  digital  signal  processing  "toolkit."  The  user 
can  select  functions  which  can  operate  on  the  signals  and  interactively  apply  them 
in  any  order.  The  user  can  also  easily  develop  new  functions  and  include  them  in 
the  "toolkit." 

The  thesis  includes  brief  discussions  about  the  library  workspaces,  a  user  man- 
ual, function  listings  with  examples  of  their  use,  and  an  application  paper.  The 
software  is  modular  and  can  be  expanded  by  adding  additional  sets  of  functions. 
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I.    INTRODUCTION 

The  Electrical  and  Computer  Engineering  Department  at  the  Naval  Postgrad- 
uate School  has  a  need  for  signal  processing  software  packages  for  use  by  students 
studying  digital  signal  processing  (DSP).  The  purpose  of  this  thesis  is  to  provide  a 
set  of  tools  in  the  form  of  APL  workspaces  on  floppy  disks  that  students  can  use  at 
home  or  at  school  on  IBM-PC-compatible  computers  for  solving  digital  signal  pro- 
cessing problems.  Workspace  is  a  collection  of  functions  and  data  that  are  placed 
together  in  order  to  do  any  job. 

The  software  supports  the  following  currently  offered  courses  in  DSP: 

1.  EC  2400-Discrete  Systems 

2.  EC  3400-Introduction  to  Digital  Signal  Processing 

3.  EC  3410-Introduction  to  Discrete  Time  Random  Processes 

4.  EC  4430- Advanced  Signal  Processing  and  Spectral  Estimation 

5.  EC  4440- Multidimensional  Digital  Signal  Processing 

Specific  APL  workspaces  were  developed  for  EC3400,  EC3410,  and  EC4440. 
The  software  uses  STSC  Corporation's  APL*PLUS  PC  and  STATGRAPHICS  sys- 
tems which  are  both  site  licensed  to  the  Naval  Postgraduate  School.  APL*PLUS 
PC  is  STSC's  version  of  the  APL  language  for  the  IBM-PC  and  compatibles.  Sim- 
ilar versions  of  the  APL  language  are  also  available  from  STSC  for  other  com- 
puters. STATGRAPHICS  (Statistical  Graphics  System)  is  a  PC  software  package 
integrating  a  variety  of  statistical  functions  with  color  graphics.  STATGRAPH- 
ICS is  written  in  APL  (specifically  APL*PLUS  PC)  and  offers  direct  access  to  the 


host  language.  Currently  STATGRAPHICS  is  available  only  for  the  IBM-PC  and 
compatibles. 

The  DSP  software  developed  in  this  thesis  is  in  the  form  of  a  digital  signal 
processing  "toolkit."  The  user  can  select  functions  to  operate  on  the  signals  and 
interactively  apply  them  in  any  order.  The  user  can  also  easily  develop  new  func- 
tions and  include  them  in  the  "toolkit."  The  tools  are  provided  on  a  set  of  floppy 
disks  that  users  can  carry  to  their  own  PC  machine.  The  sets  of  disks  contain: 

1.  The  APL*PLUS  Executive  Program 

2.  STATGRAPHICS  disks 

3.  UTILITY  workspace 

4.  EC3400  workspace 

5.  EC3410  workspace 

6.  EC4440  workspace 

The  main  products  of  the  research  are  the  disks  that  include  the  signal  pro- 
cessing library  workspace,  and  the  documentation  provided  in  this  thesis.  A  typical 
library  contains  functions  such  as:  Fast  Fourier  Transform  (FFT),  Inverse  FFT, 
Sine,  Cosine,  Sample  Autocorrelation,  Mean,  and  other  similar  functions.  In  de- 
veloping the  DSP  software  the  following  requirements  were  met: 

1.  Ability  to  input/output  from/to  a  disk  data  file  in  free  format.  This  allows 
convenient  communication  with  DOS  and  with  programs  written  in  other  lan- 
guages or  on  other  computer  systems. 

2.  Ability  to  provide  for  graphic  output  to  represent  signals,  frequency  responses, 
and  other  functions  including  2-D,  3-D,  and  contour  plotting. 

3.  Ability  to  perform  filtering,  convolution,  and  other  common  signal  processing 
operations  in  the  signal  (time/space)  domain. 


4.  Ability  to  compute  transfer  functions,  to  compute  frequency  responses,  and 
generally  to  perform  complex  arithmetic  and  other  common  frequency  domain 
operations. 

5.  Ability  to  generate  random  signals. 

6.  Ability  to  compute  statistical  characteristics  of  the  signals  including  correla- 
tion functions,  spectra,  and  prediction  errors. 

The  user  can  combine  functions  in  any  desired  order  to  solve  DSP  computer 
assignments  and  can  concentrate  on  the  signal  processing  exercises  without  the 
need  to  do  complicated  programming. 

Although  the  intent  was  to  develop  the  signal  processing  toolkit  specifically  for 
the  IBM-PC  and  compatibles,  the  system  is  portable.  STSC  has  APL  available  on 
IBM  mainframes,  Digital  Equipment  Corporation  VAX,  (under  UNIX  and  VMS), 
Apple  Macintosh,  and  others.  These  DSP  workspaces  are  portable  to  all  of  these 
systems  with  few,  if  any,  changes.  In  addition,  the  DSP  workspaces  are  easily 
carried  over  to  other  APL  systems  with  minor  changes  to  only  some  of  the  most 
basic  functions;  bits  and  pieces  of  the  software  have  existed  for  some  time  under 
VS-APL  on  the  IBM  mainframe  and  under  UNIX  4.3  bsd. 

The  remainder  of  the  thesis  is  organized  as  follows. 

Chapter  II  contains  a  brief  description  of  the  APL  language  and  the  STAT- 
GRAPHICS  software  system.  Chapter  III  deals  with  the  DSP  library  and  discusses 
the  application  of  the  four  workspaces.  The  Utility  workspace  includes  general 
functions  such  as  GETDATA  and  PUTDATA  for  the  data  input/output  to  DOS; 
the  EC3400  workspace  includes  basic  one-dimensional  DSP  functions;  the  EC3410 
workspace  includes  one-dimensional  DSP  functions  for  statistical  signal  processing; 
and  the  EC4440  workspace  includes  multidimensional  DSP  functions.  Chapter  IV 
outlines  the  conclusions,  and  the  benefits  of  using  the  entire  package  as  a  software 


toolkit.  Appendix  A  provides  a  simple  user  manual  for  the  software  package.  Ap- 
pendix B  provides  a  function  library  with  a  brief  description  of  each  function.  Ap- 
pendix C  contains  function  listings  and  samples  of  their  use.  Appendix  D  furnishes 
samples  of  computer  assignments  and  their  solutions  using  the  tools  developed  in 
this  thesis. 


II.    APL  AND  STATGRAPHICS  ENVIRONMENT 

The  signal  processing  software  is  written  in  APL  and  uses  the  APL*PLUS  PC 
interpreter  and  the  STATGRAPHICS  software  system.  Although  STATGRAPH- 
ICS provides  many  statistical  functions  in  itself,  only  the  plotting  functions  from 
STATGRAPHICS  are  used.  Following  is  a  brief  description  of  APL  and  the  STAT- 
GRAPHICS system. 

A.    APL  [1] 

Originally  APL  began  as  a  notation  to  express  mathematical  procedures.  Its 
originator,  Dr.  Kenneth  Iverson,  published  his  notation  in  1962  as  A  Program- 
ming Language  from  which  the  name  "APL"  is  derived  [2]. 

One  of  the  strongest  advantages  of  APL  is  its  ability  to  handle  complicated 
array  data-structures  without  extensive  programming.  The  language  allows  the 
definition  of  numerical  items  or  a  series  of  items  with  one  symbol.  Examples  of 
data  that  can  be  represented  are: 

•  SCALAR  -  A  single  value 

•  VECTOR  -  A  group  of  scalars 

•  MATRIX  -  A  group  of  vectors 

•  ARRAY  -  A  group  of  matrices 

The  language  makes  it  as  easy  to  perform  an  arithmetic  operation  using  a 
matrix  or  vector  argument  as  to  perform  the  operation  using  simple  numbers. 
This  characteristic  is  illustrated  in  the  following  example. 


M*"16  (define  vectors  with  6  numbers) 


12    3    4    5    6 


m-2*T.6 
N 

2    4    6    3    10    12 


n*m  (vector   multiplication) 

2    8    18    32    50    72  v  r  ' 


M«-3    3    pl9 
M 

12  3 
4  5  6 
7    8    9 


N«-3    3    p3xt.9 
M 

3       6      9 

12    15    18 
21    24    27 


(define  matrices  with  3  rows  and  3  columns) 


NxM 

3       12      27 

48     75  108  (matrix  arithmetic  operation) 

147    192    243 


Notice  that  when  two  items  are  of  the  same  dimensions,  the  arithmetic  operation  is 
performed  element  by  element.  Thus  N  x  M  is  not  the  usual  linear  algebra  product 
of  matrices.  The  latter  is  realized  by  a  composite  operation  denoted  by  +.  x. 

N+.xM 

90     108     126 

198     243    288 

306    378    450 


Scalar  arguments  are  automatically  extended  to  be  conformable  with  the  arrays. 
Thus  the  statement  M  x  2  yields 


2       4     6 
8     10     12 
14     16      18 


There  is  frequently  no  need  for  loops  or  similar  mechanisms  as  in  other  lan- 
guages such  as  Fortran.  From  the  user's  point  of  view  the  language  performs  the 
operation  on  the  elements  of  the  array  in  parallel  without  the  need  for  loops. 

Note  also  that  there  is  no  need  for  dimensions  because  APL  carries  the  dimen- 
sions of  an  item  along  with  its  definitions. 

APL  provides  a  function  to  help  the  user  to  define  or  to  determine  the  number 
of  elements  in  the  variable.  The  function  represented  by  the  Greek  character  "/>", 
is  called  "shape,"  since  it  gives  the  form  of  the  variable.  On  the  previous  page 
there  is  an  example  where  we  define  the  matrices  M  and  N  using  the  function  "p" . 
In  order  to  find  the  shape  of  any  variable,  use  the  function  "/>"  with  the  variable 
name  following  it.  The  following  is  an  example. 

pH 

3    3 

APL  has  a  large  variety  of  mathematical  functions  that  can  easily  manipulate 
the  data  and  can  handle  operations  on  a  simple  as  well  as  on  a  complicated  level. 
Thus,  it  is  possible  to  place  on  a  single  line  a  combination  of  several  functions  which 
work  together  in  a  pipeline.  All  of  this  can  be  done  interactively.  The  flexibility  in 
manipulating  array  data  easily  makes  the  language  well  suited  for  single  processing 
analysis.  The  following  are  a  few  examples: 


2    3 

3      1'1 
9      2      0 

10  3       1 

11  4       2 

4    5    <> .  +    6    "1    "3    4 
6 
7 
3 
9 

M«-3    3 

M 
13    5 
7    9    0 
3    6    4 

pi    35790354 

(outer  product)  used  to  define  separable 
2-D  function  x(nl,n2)  =  xl(nl)x2(n2) 


The  domino  operator  (calculate  M  inverse) 


BM 

"0.1313131313 
0.1414141414 
0.15 15 15 IS IS 

Oil  3  6 


"0.09090909091   0.2272727273 
0.1313131313   "0.1767676763 
"0. 09090909091   0.06060606061 


CHK 
"0.9090909091  1.04040404  1.757575753 


canj*.xc 

"0. 9090909091  1.04040404  1.757575753 


(solve  Mx  =  C  for  x) 


(check  solution  by 

multiplying  constant  vector 

by  M  inverse) 

-K  x  is  matrix  product 


APL  uses  symbols  to  represent  its  own  built-in  functions  and,  in  addition,  APL 
allows  the  user  to  define  other  functions.  The  solution  to  one  problem  consists  of 
a  collection  of  the  programmer's  functions,  each  one  of  which  can  be  applied  inter- 
actively by  the  user  or  can  be  organized  to  call  another  function.  The  advantage 
of  this  organization  is  the  ease  of  debugging  and  that  only  a  few  lines  of  APL  can 
do  volumes  of  work. 

APL  has  a  concept  called  the  workspace  which  is  a  collection  of  functions  and 
data  that  are  placed  together  in  order  to  do  a  job.  That  collection  is  stored  in  the 
computer  memory  during  an  APL  session  and  gives  the  operator  immediate 


feedback  and  accessibility,  and  the  workspace  can  be  saved  on  disk  and  loaded  back 

into  the  memory  at  any  time  for  future  use. 

B.    STATGRAPHICS  SOFTWARE  PACKAGE  [3] 

STATGRAPHICS  is  a  Statistical  Graphics  System  written  in  APL  for  scientific 
and  engineering  applications  and  for  data  analysis.  The  software  supports  dot 
matrix  printers,  monochrome/color  terminals,  and  pen  plotters.  It  runs  under  the 
APL*PLUS  PC  interpreter  version  5.0  or  later,  which  allows  the  user  to  write 
his  own  functions.  Using  STATGRAPHICS  in  this  mode  requires  familiarity  with 
APL. 

The  system  also  has  a  menu  that  contains  graphic  functions,  a  set  of  data 
analysis  and  statistical  procedures.  This  menu  is  illustrated  in  Figure  2.1.  In 
addition,  the  menu  contains  data  management  and  system  utilities  that  allow  the 
user  to  select  the  system  environment,  control  color,  plot  dimensions,  and  so  on. 
The  most  often  used  screens  for  DSP  applications  are: 

1.  Plotting  functions  (Figure  2.2) 

2.  Time  Series  Analysis  (Figure  2.3) 

3.  Experimental  Design  (Figure  2.4) 

STATGRAPHICS  has  options  for  modifying  graphs  and  for  on-line  help. 

To  use  STATGRAPHICS  with  APL  the  user  needs  to  start  APL*PLUS  PC  as 
usual  and  then  load  the  STATGRAPHICS  workspace  from  the  start-up  disk.  By 
going  back  and  forth  from  the  APL  environment  to  STATGRAPHICS,  the  user 
can  process  data  and  plot  the  results  easily.  This  makes  the  combination  of  APL 
and  STATGRAPHICS  powerful  and  very  convenient  to  use. 


STATGRAPHICS  Statistical  Graphics  Systami 


DATA  MANAGEMENT  AND  SYSTEM  UTILITIES 

A.  Data  Management 

B.  System  Environment 

C.  Report  Writer  and  Graphics  Replay 

D.  Plotter  Interface 

PLOTTING  AND  DESCRIPTIVE  STATISTICS 
S.  Plotting  Functions 

F.  Descriptive  Methods 

G.  Estimation  and  Testing 
H.  Distribution  Functions 

I.  Exploratory  Data  Analysis 

ANOVA  AND  REGRESSION  ANALYSIS 
J.  Analysis  of  Variance 
X.  Regression  Analysis 


TIJCE  SERIES  PROCEDURES 
L.  Forecasting 
M.  Quality  Control 
tf.  Smoothing 
0.  Time  Series  Analysis 

ADVANCED  PROCEDURES 

P.  Categorical  Data  Analysis 

Q.  Multivariate  Methods 

P..  Nonparametric  Methods 

S.  Sampling 

T.  Experimental  Design 

MATHEMATICAL  AND  USER  PROCEDURES 
U.  Mathematical  Functions 
V.  Supplementary  Operations 


Use  cursor  Key a  tu  highlight  desired  section.   Then  press  ENTER. 
IHelp   2Edit   3Savscr  4Prtscr  5       SGo     7Vars   3Cmd    9Review  lOQuit 
INPUT     11/13/37   13:03  STATGRAPHICS  Vers.  2.1 


Figure  2.1.    STATGRAPHICS  Main  Menu. 


PLOTTING  FUNCTIONS 


1.  2-Y  Line  and  Scatterplots 

2.  Multiple  X-Y  Plot3 

3.  X-Y-Z  Line  and  Scattarplota 
*.  Multiple  X-Y-Z  Plots 

5.  Barcharts 

S.  Piecharts 

7.  Component  Line  darts 


Use  cursor  keys  to  highlight  desired  procedure.   Then  press  ENTER. 
IHelp   2Edit   3Savscr  4Prtscr  5       SGo     7Var3   3Cmd    9Review  lOQuit 
INPUT    11/13/37   13:10  STATGRAPHICS  Vers.  2.1 


Figure  2.2.    The  Plotting  Function  Menu. 


10 


TIME  SERIES  ANALYSIS 


1.  Horizontal  Time  Sequence  Plot 

2.  Vertical  Tine  Sequence  Plot 

3.  Seasonal  Subseries  Plot 

4.  Autocorrelation  Function 

5.  Partial  Autocorrelation  Function 

6.  Cross -Correlation  Function 

7.  Simple  or  Seasonal  Differencing 
3.  Mean  or  Trend  Removal 

9.  Box-Cox  Transformation 

10.  Periodogram 

11.  Integrated  Periodogram 

12.  Data  Tapering 

13.  Plotting  vs.  Fourier  Frequencies 

14.  Sox- Jenkins  ARIMA  Modeling 

15.  Cross-Correlation  Matrix  Plot 


Use  cursor  keys  to  highlight  desired  procedure. 
lHelp   2Edit   3Savscr  4Prtscr  5       6Go     7Vars 
INPUT     11/18/87   13:13  STATGRAPHICS  Vers.  2.1 


Then  press  ENTER. 
8Cmd    9Review  lOQuit 


Figure  2.3.    The  Time  Series  Analysis  Menu. 


EXPERIMENTAL  DESIGN 


1.  Full  and  Fractional  Factorials 

2.  Central  Composite  Designs 

3.  Alias  Structure 

4.  Response  Surface  Plotting 


Use  cursor  keys  to  highlight  desired  procedure. 
lHelD   2Edit   3Savscr  4Prtscr  5       6Go     TVars 


Then  press  ENTER. 
8Cmd    9Review  lOOuitr 


Figure  2.4.    The  Experimental  Design  Menu. 
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III.    THE  DSP  WORKSPACES  AND  APPLICATIONS 

The  following  is  the  principle  part  of  the  thesis  research.  This  chapter  describes 
the  workspaces  that  were  developed  by  the  author,  and  the  way  they  can  be  applied. 

A.    UTILITY  WORKSPACE 

This  workspace  contains  functions  of  general  use  that  can  be  applied  in  combi- 
nation with  functions  from  other  workspaces  in  this  package.  Appendix  B  contains 
the  list  of  workspace  functions  and  a  brief  discussion  of  each  one. 

Two  of  the  most  useful  functions  in  the  utilities  are  described  and  demon- 
strated below.  The  functions  GETDATA  and  PUTDATA  allow  the  user  to  transfer 
data  from  the  APL  environment  to  the  DOS  environment  and  vice  versa.  These 
functions  use  file  system  commands  that  are  specific  to  STSC's  implementation  on 
the  IBM-  PC.  As  such,  the  exact  APL  code  for  GETDATA  and  PUTDATA  is  not 
always  portable  to  APL  on  other  systems.  However,  versions  of  GETDATA  and 
PUTDATA  can  be  written  for  any  system.  Versions  are  currently  available  for  the 
IBM  mainframe  under  VS/APL  and  for  a  version  of  APL  under  UNIX  4.3  bsd  on 
the  VAX.  To  the  user,  GETDATA  and  PUTDATA  provide  a  common  interface  for 
communication  with  other  programs. 

The  following  is  an  example  showing  the  transfer  of  data  using  the  GET- 
DATA and  PUTDATA  functions.  Figure  3.1  describes  two  sets  of  data  in  an  APL 
workspace.  In  order  to  transfer  those  sets  to  a  file  in  the  DOS  environment  the 
user  needs  to  type  the  following  commands: 
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X  PUTDATA    'X.DAT' 
Y  PUTDATA    'Y.DAT' 
X  is  5  x  5  and  Fisa2x4x4  array 
x 

0.925  1.35  2.775  3.7  4.525 

5.55  6.475  7.4  3.325  9.25 

10.175  11.1  12.025  12.95  13.375 

14.3  15.725  1S.S5  17.575  13.5 

19.425  20.35  21.275  22.2  23.125 

Y 

1ES0  2ES0  3E50  4ES0 

5ES0  5250  7ES0  3ES0 

9ES0  1ES1  "1E44  "2E44 

"3E44  "4E44  "5E44  "SE44 

"7E44  "3E44  "9E44  "1E45 

"1.1E45  "1.2E4S  "1.3E45  "1.4E45 

"1.5E45  1.SE45  "1.7E4S  "1.3E4S 

"1.9E45  ~2E45  "2.1E45  "2.2E45 

Figure  3.1.    X  and  Y  Data  Sets. 

The  left  argument  to  PUTDATA  is  the  variable  name;  the  right  argument  is  the 
file  name  or  path.  Note  that  after  transferring  any  data  set  to  DOS  the  data  is  not 
erased  and  remains  available  in  the  workspace. 

When  terminating  APL,  the  user  is  placed  in  the  DOS  environment.  Figure 
3.2  shows  that  the  data  set  transfer  was  successfully  completed.  Note  that  the  first 
line  of  the  data  set  gives  the  dimension  or  "shape"  of  the  data. 

By  going  back  to  APL  environment  and  typing  the  commands: 

X<-GETDATA    'X.DAT' 
Y<-GETDATA    'Y.DAT' 
the  data  is  transferred  correctly  from  the  DOS  environment  to  the  APL  workspace. 
Figure  3.3  shows  that  the  data  has  the  same  original  shape.    The  data  are  read 
into  the  workspace  and  assigned  to  the  variables  X  and  Y .  The  data  still  remain 
available  as  a  file  in  the  DOS  environment. 
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)OFF 


C:\APL>TYPE   X.DAT 


(terminate  APL) 


(type  of  contents  of  files) 


5 

5 

0.925 

1.35 

2.775 

3 

.7 

4.625 

5.55 

6.475 

7 

.4 

3.325 

9.25 

10.175 

11 

.1 

12.025 

12.95 

13.375 

14 

.3 

15.725 

16.65 

17.575 

13 

.5 

19.425 

20.35 

21.275 

22 

.2 

23.125 

:\APL>TYPE  Y.DAT 

2 

4 

4 

1E60 

2E60 

3E60 

4E60 

5ES0 

6E60 

7E60 

3E60 

9E60 

1E61 

-1E44 

-2E44 

-3E44 

-4E44 

-5E44 

-6E44 

-7E44 

-8E44 

-9E44 

-1E45 

-1.1E45    -1.2E45    -1.3E45    -1.4E45 
-1.5E45    -1.SE45    -1.7E45    -1.3E45 
-1.9E45    -2E45    -2.1E45    -2.2E45 
C:\APL> 


Figure  3.2.    X  And  Y  Data  Sets  In  DOS  Environment. 


pX 

5  5 

X 

0.925   1.35    2. 775 

3.7     4.625 

5.55    6.475   7.4 

3.325   9.25 

10.175  11.1    12.025 

12.95   13.375 

14.3    15.725  16.55 

17.575  13.5 

19.425  20.35   21.275 

22.2    23.125 

P* 


2    4    4 


1E60 

2E60 

3E60 

4E60 

5E60 

6E60 

7E60 

3E60 

9E60 

1E61 

"1E44 

"2E44 

"3E44 

"4E44 

"5E44 

"6E44 

"7E44 

"3E44 

"9E44 

"1E45 

"1  1E45 

"1.2E45 

"1.3E45 

"1.4E45 

"1.5E45 

"1.6E45 

"1.7E45 

"1.3E45 

"1.9E45 

"2E45 

"2.1E45 

"2.2E45 

Figure  3.3.    Getting  X  And  Y  Into  APL  Workspace. 

B.    EC3400  WORKSPACE 

EC  3400  is  a  course  at  the  Naval  Postgraduate  School  that  introduces  first 
principles  of  digital  signal  processing.    The  topics  covered  include:    the  Discrete 
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Fourier  Transform  and  the  FFT  algorithm,  flow-graph  and  matrix  representation 
of  filters,  ideal  filters  and  approximation,  and  design  of  recursive  and  nonrecursive 
digital  filters,  fast  convolution  and  correlation. 

The  following  workspace  enables  the  student  to  solve  computer  assignments  on 
the  above  subjects  in  a  simple  way.  This  workspace  contains  functions  such  as:  lin- 
ear convolution,  circular  convolution,  FFT,  inverse  FFT,  and  complex  arithmetic. 
Some  of  those  functions  will  be  described  later.  In  addition,  the  EC3400  workspace 
contains  functions  which  allows  the  student  to  design  one- dimensional  filters.  The 
methods  used  are  the  Fourier  method  and  the  frequency  sampling  method  and 
windowing  of  the  filters  coefficients.  The  full  description  of  the  workspace  is  given 
in  Appendix  B.  The  following  is  an  example  of  how  to  use  this  workspace. 


EC  3400  Computer  Assignment  [4 


A  highpass  filter  is  to  have  an  analog  cutoff  frequency  of  6kHz.  The  filter  is 
to  be  implemented  by  a  digital  filter  having  a  sampling  frequency  of  40kHz.  Plot 
the  frequency  response  of  the  filter  with  and  without  a  Hamming  window  if  51 
coefficients  are  used. 


Solution 

The  solution  given  here  uses  functions  provided  in  the  EC3400  workspace.  In 
an  actual  classroom  assignment,  students  may  be  asked  to  generate  their  own  func- 
tions or  to  solve  the  problem  in  their  own  way.  By  using  the  function  DIGFREQ, 
the  digital  cutoff  frequency  is  calculated.  (Figure  3.4).  This  function,  like  several 
others  is  included  for  convenience  only;  the  calculation  in  converting  radian  to  dig- 
ital frequency  can  be  done  easily  in  the  workspace  by  using  the  built-in  arithmetic 
operations.     Figure  3.5  describes  the  ideal  frequency  response  to  be  generated. 
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By  using  the  function  HPCOEFF,  51  causal  impulse  response  coefficients  are  cal- 
culated. (See  Figures  3.6  and  3.7).  The  function  HPCOEFF  uses  the  following 
equation: 

hLP{n)  =  k/ir(n  -  I)  sin([n  -  I]0C),      n  =  0, 1, 2, . . . ,  21 


TC«-40000  DIGFBEQ  3000 
TC 
1.2i»663  7061 


TC+PI 


0.4 


Figure  3.4.    Calculating  a  Digital  Cutoff  Frequency. 
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Figure  3.5.    The  Ideal  Frequency  Response. 
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hHP«-51  HPCOEFF  TC 

hHP 
6.363316738E"17  "0.01261377881  "3.13 
2.339861191E~17  "0.01593319429 
0.01892066822  2. 339861191E"17 
0.01559148806  0.02752097195  2 
0.02672826525  0.05045511524  2 
"0.06236595225  0.09354892838  0 
0.09354892838  "0.06236595225  " 
0.05045511524  0.02672826525  "0 
0.02752097195  0.01559148806  "0 
2.339861-191E~17  0.01392066322 
"0.01593319429  2.339861191E-17 
~3.134689424E"3  "0.01261377881 


4689424E~3  3. 504448034E"3  0.01441574 

"0.01039432538  0.01100575628 
"0.02162362082  "0.01439214283 
339861191E"17  "0.0336367435  "0.02338 
339861191E"17  "0.07563267236 

3027306915  "0.6  0.3027306915 
0.07563267286  2. 339861191E~17 
.02338723209  "0.0336367435  2.3393611 
.01439214283  "0.02162362082 
0.01100575623  "0.01039432538 

0.01441574721  3.504448034E"3 

6.363316738E"17 


Figure  3.6.    Calculating  the  Impulse  Response  Coefficients. 
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Figure  3.7.    The  Impulse  Response  Coefficients. 

The  function  HAMMING  generates  "J"  samples  of  a  shifted  Hamming  window 
using  the  following  equation: 


W(n)  =  0.54  +  0.46  cos 


2ttti 
I-l 


,0  <n<  I-l 


This  is  multiplied  by  the  filter  coefficients  to  obtain  the  windowed  coefficients.  (See 
Figures  3.8  and  3.9). 

By  using  the  function  FREQRES,  the  frequency  response  with  and  without  a 
window  is  generated.  This  is  done  by  calculating  the  DFT  of  the  filter  coefficients 
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W<-HAMMING  51 


W 
0.08087246878  0.0878 
0.2197783849  0 
0.5258342731  0 
0.7940878876  0 
0.9689372255  0 
0.9689372255  0 
0.7940878876  0 
0.5258342731  0 
0.2197783849  0 
0.08783237415 


3237415  0.1016466798  0.122 


2627880673 
.5824434454 
.8392844181 
.9861065906 
.9452656094 
.7450396437 
.4694398388 
.1816229357 


.31  0.3606984 
.6384092183  0 
,8799441019 
,9965134344 
,9154505797 
6928832078 
4141150246 
1489001176 


105975  0.148 
983  0.414115 
.6928832078 
.9154505797 

0.996513434 
.8799441019 
.6384092183 
.3606984983 
.122105975  0 


9001176 
0246  0. 
0.74503 
94526 
0.986 
33928 
58244 
31  0. 


.101646 


0.1816 
4694398 
96437 
56094 
1065906 
44181 
34454 
2627880 
6798 


0.08087246878 


hHPW«-hHPxW 


Figure  3.8.    Calculating  the  Hamming  Window  Samples. 
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Figure  3.9.    The  Hamming  Window. 

and  finding  the  magnitude  of  the  result  as  shown  in  a  section  of  the  FREQRES 
function  (Figure  3.10).  The  following  describes  the  APL  commands  which  generate 
the  frequency  responses  and  Figures  3.11  and  3.12  show  the  corresponding  results. 

HHP*-128  FREqRES  hHP 
HHPW<-128  FREQRES  hHPW 
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CO]  HW«-M  FREQRES  h;DIO;N:H 

CI]  « 

C2]  A   GENERATING  FREQUENCY  RESPONSE  MAGNITUDE  AND  PHASE  FROM  TIME  DOMAIN 

C3]  fl   COEFFICIENTS 

[4]  A 

C5i   nio«-o 

C6]    N«-CM-ph)pO 
I2J h^-h.N 


C3]    H<-FFT  h 

C9]    HW«-CM-5-2)tCXMAGN  H) 

[10]   PW<-CM*2)TCXPHAS  g] 


C11J 

C12] 
C13] 
C14] 
C153 
[16] 
C17] 
[13] 


•+Q 


THE  FUNCTIONS  FFT . XMAGN . XPHASE  ARE  USED. 

IN  ORDER  TO  RECEIVE  A  BETTER  RESOLUTION  THE  COEFFICIENTS  ARE  ZERO- 

PADDET  TO  'M'  SAMPLES. 

M<-NUMBER  OF  SAMPLES  INCLUDING  ZERO  PADDING  CUSE  ONLY  RADIX  TWO  NO 

h<-THE  FILTER  COEFFICIENTS  [VECTOR: 

HW: FREQUENCY  RESPONSE  MAGNITUDE  CM*  2  SAMPLES) 

PW:  FREQUENCY  RESPONSE  PHASE  CM*  2  SAMPLES)  [RAD.] 


Figure  3.10.    Using  the  Function  FFT. 
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Figure  3.11.    The  Frequency  Response  Using  51  Coefficients. 

C.    EC3410  WORKSPACE 

EC  T5410  is  a  course  that  gives  an  introduction  to  discrete-time  random  process. 
The  topics  covered  are:  description  of  discrete-time  random  signals  and  random 
vectors,  linear  transformations,  sampling  of  continuous-time  random  signals,  esti- 
mation, and  spectral  analysis.  The  EC3410  workspace  was  developed  in  order  to 
enable  the  student  to  solve  computer  assignments  in  digital  signal  processing  and 
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Figure  3.12.    The  Frequency  Response  Using  Hamming  Window. 

learn  to  write  algorithms  in  an  efficient  way  without  doing  excessive  programming. 
This  workspace  contains  functions  that  compute  the  mean,  circular  convolution, 
linear  convolution,  fast  Fourier  transform,  power  spectrum  and  other  quantities. 
A  full  description  of  the  workspace  functions  is  given  in  Appendix  B. 

By  using  the  STATGRAPHICS  software  package  with  the  EC3410  workspace, 
a  graphic  capability  is  achieved.  The  following  is  a  computer  assignment  in  statis- 
tical dsp  and  its  solution  using  the  EC3410  workspace  and  the  STATGRAPHICS 
software. 

Note  that  in  the  actual  assignment,  students  are  required  to  write  some  of  the 
functions  that  are  used  below.  However  developing  the  functions  in  APL  using  the 
tools  provided  is  not  difficult. 
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EC  3410  COMPUTER  ASSIGNMENT  [5] 

1.  Three  data  sets  SL1.DAT,  SL2.DAT,  and  SL3.DAT  are  given  on  a  disk.  The 
mean  of  the  signal  data  SL1.DAT  is  0.75.  Subtract  this  constant  value  from 
the  signal  to  generate  a  new  zero-mean  signal  SL1'  that  you  will  use  in  this 
assignment.  The  other  data  sets  (SL2,  SL3)  have  zero  mean.  Do  the  remaining 
parts  for  each  of  the  data  sets  SL1',  SL2,  and  SL3. 

2.  Using  the  correlation  method,  generate  a  3  by  3  Toeplitz  correlation  matrix 
for  the  signal  data.  Print  this  matrix. 

3.  Solve  a  set  of  Normal  equations  involving  the  correlation  matrix  that  was  just 
generated  to  obtain  the  2nd  order  linear  predictive  filter  parameters  and  the 
prediction  error  variance. 

4.  Apply  the  filter  to  the  original  data  set  and  generate  the  prediction  error 
signal.  Plot  this  signal  and  compute  its  variance.  Does  it  compare  well  with 
the  theoretical  prediction  error  variance  you  obtained  by  solving  the  normal 
equations? 


Solution 

1.  By  using  the  function  MEAN,  the  SL1  data  set  mean  is  computed  and  is 
subtracted  from  the  original  SL1.  The  new  data  set  SL1'  has  a  zero-mean. 
(Figure  3.13) 


M«-MEAN 

SL1 

M 

2. 

065552043 

SL1«-SL1-M 

M«-MEAN 

SL1 

M 

6. 

982261991E" 

■17 

Figure  3.13.    Using  the  Function  MEAN. 
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Rl<-3    SACF   SL1 

Rl 
3.575461794    3.055195789    7.537100879 

Kl<-COV   Rl 

Kl 

8.575461794  3.055195739  7.537100379 

3.055195789  3.575461794  3.055195739 

7.537100879  3.055195739  3.575461794 


Figure  3.14.    Using  the  Functions  SACF  to  Evaluate  the  Covariance  for 
SL1.DAT. 


R2«-3    SACF   SL2 

R2 

1.214143223    0.589008547    0.2762704521 

K2*-C0V  R2 

K2 

1.214143223  0.589008547  0.2762704521 

0.589008547  1.214143223  0.589008547 

0.2762704521  0.589008547  1.214143223 


R3«-3  SACF  SL3 

R3 

9.410119655  "3.931635734  3.467433519 

K3*-C0V  R3 

K3 

9.410119655  "3.931635784  8.467433519 

"3.931635734  9.410119655  "3.931635784 

8.467433519  "3.931635784  9.410119655 


Figure  3.15.    Using  the  Functions  SACF  and  COV  for  SL2  and  SL3. 

2.  By  using  the  function  SACF  (sample  autocorrelation  function),  the  first  three 
correlation  function  lags  are  calculated  and  then  by  using  the  function  COV  a 
Toeplitz  matrix  is  generated  from  the  correlation  function  data.  The  process 
described  in  the  above  paragraph  is  repeated  on  the  SL2  and  SL3  data  sets, 
as  seen  in  Figures  3.14  and  3.15. 

3.  The  function  LPFP  (linear  prediction  filter  parameters)  solves  the  Normal 
equations  and  it  calculates  the  2nd  order  filter  coefficients  and  the  prediction 
error  as  demonstrated  in  Figure  3.16.   (This  could  also  be  done  without  the 
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al«-2   LPFP  SL1 


a2<-2   LPFP  SL2 


al 
1  "0.96669S1989  0.02913175079 

PERR 
1.008111S43 


Al 

1  "0.4900717314  0.01020142115 

PERR 
0.9283051353 


a3«-2  LPFP  SL3 


a3 


1  0.9593733731  0.01077400633 

PERR 
0.9325296297 

Figure  3.16.    Solving  the  Normal  Equations  Using  the  Function  LPFP. 

LPFP  function  by  setting  up  simple  linear  equations  and  solving  them  as 
described  in  Chapter  II.) 

4.  Using  the  LCV  (linear  convolution)  function,  the  three  data  sets  are  convolved 
with  the  2nd  order  filter  parameters  that  were  found  in  the  previous  step  to 
"  provide  the  prediction  error  sequence.  The  prediction  errors  of  the  signals 
are  determined  by  computing  the  initial  value  of  the  autocorrelation  function 
of  the  prediction  error  sequence  as  seen  in  Figure  3.17.  As  can  be  seen,  the 
results  are  very  similar  to  the  prediction  error  values  derived  by  the  analysis 
of  Figure  3.16.  Figure  3.18  is  a  plot  of  the  prediction  sequence  for  each  data 
set. 

Additional  examples  of  using  the  EC3410  workspace  are  provided  in  Appendix 
D. 


D.    EC4440  WORKSPACE 

EC  4440  is  a  course  at  the  Naval  Postgraduate  School  in  multidimensional 
digital  signal  processing.  The  course  deals  with  the  analysis  of  signals  that  are 
functions  of  the  two  or  more  independent  variables.  The  course  develops  both 
time/space  and  frequency  domains  approaches  and  highlights  the  subjects  that 
are  different  from  one-dimensional  signal  processing.   Some  of  the  topics  covered 
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al«-2  LPFP  SL1 


1  "0.9666951989  0-02913175079 

PERR 
1.008111643 


a2<-2  LPFP  SL2 


a2 

1  "0.4900717314  0.01020142115 

PEER 
0.9233051353 


a3*-2  LPFP  SL3 


a3 
1  0.9593733731  0.01077400633 

PERR 
0.9325296297 


Figure  3.17.    Calculating  The  Prediction  Error  Using  The  Function  PRER. 

are:  two-dimensional  circular  and  linear  convolution,  difference  equations,  recur- 
sively computable  systems,  Fourier  analysis,  and  methods  of  2-P  filter  design.  The 
EC4440  workspace  enables  the  student  to  solve  computer  assignments  in  multidi- 
mensional signal  processing  by  using  functions  such  as:  2-D  convolution,  2-D  Fast 
Fourier  Transform  (FFT)  and  inverse  2-D  FFT,  and  filter  design  using  the  win- 
dowing method  or  the  McClellan  transformation.  The  plots  are  provided  by  using, 
as  before,  the  STATGRAPHICS  software  package.  More  detailed  lists  and  descrip- 
tions of  the  workspace  functions  are  provided  in  Appendices  B  (DSP  library)  and 
C  (functions  listings  and  samples  of  runs). 
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PREDICTION  ERROR  OF  SL1 


4.  4- 


E 
R 
R 
O 
R 


lOO   200   300   -4-00   SOO   ©OO 

TIME  SfiMPELS 


PREDICTION  ERROR  OF  SL2 


E 
R 
R 
O 
R 


O     lOO    200   300   400    SOO    £00 
TIME  SAMPLES 


PREDICTION  ERROR  OF  SL3 


E 
R 
R 
O 
R 


lOO    200   300   -4-00    500    &00 

TI ME  SAMPLES 


Figure  3.18.    The  Prediction  Errors  For  SLl,  SL2,  And  SL3. 
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The  following  is  an  example  of  how  to  design  a  lowpass  filter  using  the  workspace 
functions.     The  design  here  is  done  by  the  McClellan  transformation  method. 
Again,  as  in  other  examples,  the  students  are  asked  to  write  some  of  these  functions. 


COMPUTER  ASSIGNMENT  ON  McCLELLAN  TRANSFORMATION  [6] 
Consider  the  prototype  filter: 


BM 


—  7T  —  A-K  Alt  IT  Ui 


1.  Using  a  32  point  FFT,  find  the  impulse  response  hp(n)  and  the  coefficients 
a(n)  in  the  representation 


N 


Hp(u)  =  ^T  a(n)Tn[cosoj]         TV  =  15 

2.  Now  let  F(ll>i,uj2)  =  \  (—1  4-cosu;i  +  cosu;2  +  coswi  •  cosu^) 

Compute  H((jJi,ui2)  =  Yln=o  a{n)r^n[F{ijJiiljJ2)\  a^d  generate  both  3-D  and 
contour  plots  of  this  frequency  response. 

Note:    You  can  take  advantage  of  the  recursion 

Tn(x)  =  2xTn_1(x)-Tn_2(x) 

and  the  built-in  recursive  nature  of  APL  functions  to  generate  the  Chebyshev 
polynomials. 
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Solution 

The  function  PROFILT  is  used  to  specify  a  32-sample  ideal  lowpass  prototype 
filter;  then  by  using  the  function  COEFF  the  impulse  response  coefficients  are 
calculated.  As  shown  in  the  following  steps  the  COEFF  function  uses  the  IFFT 
function  to  calculate  the  impulse  response,  h(n),  of  the  prototype  filter. 


7C0EFFCCn? 

CO]  a«-coEFF  Hp:aiO:SlZZ:HH:H 

CI]  a 

C21  A      COMPUTING  THE   COEFFICIENTS    FROM  THE   PROTOTYPE    FILTER    C USING    IFFT] 

C3]  A 

C4]  QIO«-0 

C5]  SIZE*-(pHp)+2 

Z61  a<-SIZEoO 


C7j         H<-IFFT    Hp 


a(n)  = 


The  filter  coefficients  are  defined  as: 

A(0),       n  =  0 

2h(n),     n  >  0 

and  are  found  by  performing  the  following  steps  in  the  function  COEFF. 


C33         HH«-(1.  SIZED  TH 
191        HH«-SIZEoHH 


C10J      aCi.SIZE]«-2xHHCl.SIZE] 
Cll]      aCQ]<-aCQ]*2 


C12J  -*Q 

C13]  A 

C14]  A   Hp«-THE  PROTOTYPE  SAMPLES 

CIS]  A   a: THE  IMPULSE  RESPONSE  COEFFICIENTS 

CIS]  A   THE  IFFT  FUNCTION  HAS  BEEN  USED. 

C17]  A 

C13]  A   Y.KATZIR.I.A.F.  .  OCTOBER  1987 

C19]  A 
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HP«-PROTFILT  32 


a«-COEFF  HP 


HP 
11111110000000000000000000111111 


0.40625  0.610186859  0.1779849416  "0.1365886089  "0.1508883476  0.0129955742 
0.110335429  0.0464416897  "0.0625  "0.07130573562  0.01466457095 
0.07052675554  0.02588834765  "0.0504871114  "0.05298494156  0.01823057754 

Figure  3.19.    Calculating  the  Impulse  Response  Coefficients. 
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Figure  3.20.    1-D  Prototype  Filter. 

The  above  calculations  are  demonstrated  in  Figure  3.19.  (Notice  that  the  prototype 
filter  response  is  defined  on  the  interval  0  to  2x  rather  than  —  x  to  4-7T.)  Figures 
3.20  and  3.21  are  plots  of  the  prototype  filter  and  the  prototype  filter  coefficients. 
The  ideal  filter  would  have  a  step  discontinuance  at  the  points  u  =  0.47T 
and  1.6tt  but  this  appears  in  the  plot  as  slanted  fines  because  the  plot  uses  finite 
numbers  of  samples. 
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PS070TYPI  FILTS2  , COEFFICIENTS  h<n) 
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Figure  3.21.    The  Prototype  Filter  Coefficients. 

The  function  TRANSFNC  generates  the  transformation  function  by  imple- 
menting the  following  equation: 

F(wl,  w2)  =  0.5(— 1  4-  cos  wl  4-  cos  w2  4-  cos  tul  •  cos  u/2). 

The  surface  plot  capability  of  STATGRAPHICS  is  used  to  form  a  3-D  plot  and  a 
contour  plot  of  the  transformation  function  as  seen  in  Figure  3.22. 

By  using  the  function  MCCLEL  the  frequency  response  of  the  desired  filter  is 
generated.  The  function  generates  the  nth  Chebyshev  polynomial,  Tn[F(wl,  w2)] 
by  using  the  recursion  Tn(x)  =  2(x)Tn-i(x)  —  T„_2(x).  By  using  the  following 
summation: 

N 

£  a(n)Tn  [F(wl,  w2)] 


n=0 
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FW-TBANS7NC  32 
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Figure  3.22.    The  Transformation  Function. 
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the  frequency  response  H(wl,w2)  is  achieved.    The  above  steps  are  shown  in  a 
section  of  the  MCCLEL  function. 


VMCCLELCD]?  

CO]  Hww«-a  MCCLEL  F;niO;N:INTER:NUMBER:INTERl:INTERO 

CI]  a 

C2]  A   DESIGNING  2-D  FILTER  USING  MEcCLELLAN  TRANSFORMATION 

C3]  A 

C4]  nio*o   

C5]  INTERO*-INTERl«-INTER«-Hwv+'CpF]pO 

C6]  NUMBER«-oa        


L7J    H*-Q 

C8]  LOOP:-»(  CN=03  .N=13 /A.B 

C9]    INTER«-(2xFxINTER13-INTER0 

CIO]  -*C 

Cll]  A:INTER«-CpF]pl 

C12]   INTER1*- INTER 

C13]   -»C 

C14]  B:INTER«-F 

C15]  C:  Q«-'  CALCULATED  CHEB.  POLYNOM  OF  ORDER:  '.  *N 

C16]   Hww«-HwvKaCN]xINTER] 

C17]   INTER0«-INTER1 

C13]   INTER1«-INTER 

C19]   -»(NUMBER>N<-N-*-l)/LOOP 


L20J 
C21] 
C22] 
C23] 
C24] 
C25] 
C26] 
C27] 
C28] 


-»0 


F«-THE  TRANSFORMATION  MATRIX  (CALCULATED  USING  TRANSFORMATION 

FUNCTION] 

a«-THE  PROTOTYPE  IMPULSE  RESPONSE  COEFFICIENTS 

HwwiTHE  2-D  FREQUENCY  RESPONSE. 

Y.KATZTR.I.A.F. .  OCTOBER  1987 


Figure  3.23  shows  the  frequency  response  plot  of  the  lowpass  filter.  In  order 
to  apply  a  status  report  during  calculation,  a  status  printout  was  added  to  the 
function  as  seen  in  Figure  3.24. 

E.    DISCUSSION 

The  workspaces  and  functions  developed  in  this  thesis  provide  a  complete 
set  of  tools  for  solving  the  problems  currently  assigned  to  students  in  courses  EC 
3400,  EC  3410,  and  EC  4400.  They  also  provide  capabilities  to  solve  many  other 
problems  in  DSP. 
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HUiiK-a  MCCLIL  FU 


40 


Contour  plot  CHCcLZLLAN  Transformation) 


'      id 


Figure  3.23.    The  Frequency  Response  of  the  Lowpass  Filter. 
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Figure  3.24.    Status  Printout  During  Calculation  of  Chebyshev  Polyno- 
mial. 

Some  of  the  functions  developed  and  included  in  the  workspace  are  finite 
specialized  and  were  specifically  coded  as  examples  of  solving  the  particular  DSP 
homework  assignments.  Others  are  provided  as  tools  that  can  be  applied  to  many 
different  problems.  As  already  mentioned,  it  is  very  easy  to  expand  this  set  of  tools 
to  include  others. 

In  using  the  workspaces  and  functions,  the  instructor  will  usually  provide 
only  a  subset  of  the  functions  developed  here  to  students  enrolled  in  the  course. 
Students  would  then  use  the  functions  provided,  as  well  as  the  general  features  of 
APL  to  solve  the  assigned  problems.  Since  the  workspaces  are  extremely  modular, 
an  instructor  will  find  it  easy  to  select  only  those  functions  that  he  wants  students 
to  have  and  provide  others  at  a  later  time. 
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IV.    SUMMARY  AND  CONCLUSION 

The  purpose  of  this  thesis  was  to  develop  a  convenient  set  of  software  tools  for 
students  and  instructors  in  order  to  solve  computer  assignments  in  digital  signal 
processing  on  personal  computers  using  APL. 

The  main  requirement  was  achieved  by  developing  APL  workspaces  with  sets  of 
functions  that  could  be  applied  to  perform  DSP  operations.  Computer  assignments 
in  EC  3400.  EC  3410,  and  EC  4440  are  easily  solved  by  using  this  software  package. 
The  user  needs  to  have  basic  knowledge  in  APL  and  needs  to  know  how  to  use  part 
of  the  STATGRAPHICS  package  in  the  APL  environment.  However,  this  basic 
knowledge  is  easily  acquired. 

The  structure  of  the  APL  software  is  very  flexible  and  convenient  to  use.  By 
combining  a  few  functions,  the  user  can  achieve  complicated  calculations  with  little 
programming  effort.  Any  workspace  can  be  expanded  by  the  student  or  professor 
by  adding  new  functions. 

The  results  of  this  thesis  demonstrate  that  the  APL  is  easy  and  convenient 
to  use,  especially  in  solving  digital  signal  processing  problems.  This  thesis  will 
hopefully  make  it  easier  for  new  students  to  start  using  the  language. 
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APPENDIX  A 
USER  MANUAL 

A.  SYSTEM  REQUIREMENTS 

The  software  package  contains  the  following  disks: 

1.  APL*PLUS  Executive  Program  version  5.1  or  later  by  STSC 

2.  STATGRAPHICS  -  Statistical  Graphics  System,  version  2.1  or  later  by  STSC 

3.  Utility  workspace,  and  EC3410,  EC4440,  and  EC3400  workspaces 

STATGRAPHICS  requires  the  following  hardware  and  system  software: 

1.  IBM  PC,  PC-XT,  PC- AT  or  compatible 

2.  512K  RAM 

3.  Keyboard  with  APL  characters 

4.  Two  double-sided,  double-density  disk  drives  or  one  floppy  disk  and  hard  disk. 

5.  Graphics  adapter  card 

6.  DOS,  version  2.0  or  later 

7.  Dot  matrix  graphics  printer  (optional) 

8.  Math  coprocessor  (optional  but  recommended) 

The  software  package  can  be  operated  from  floppy  disks  or  from  the  hard  disk. 

B.  PREREQUISITES  FOR  USER 

In  order  to  use  the  package  efficiently,  the  user  should  have  a  basic  knowledge 
of  APL. 

In  addition,  it  is  strongly  recommended  that  a  copy  of  the  STATGRAPHICS 
manual  [3]  be  available  to  the  user  since  this  Appendix  refers  the  user  to  Ref.  3 
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for  more  explanations.  The  user  should  not  have  loaded  any  special  DOS  resident 
programs  since  this  may  interfere  with  the  proper  execution  of  STATGRAPHICS 
and  cause  the  system  to  "hang  up." 

C.    HOW  TO  GET  INTO  THE  DSP  LIBRARY 

Follow  the  instructions  below  to  start  the  session: 

1.  Start  the  computer  using  DOS. 

2.  If  you  are  using  a  dual  floppy  disk  drive,  place  the  APL*PLUS  executive  disk 
in  the  default  disk  drive.  If  you  are  using  a  hard  disk  system,  make  sure  that 
the  directory  with  the  APL  and  STATGRAPHICS  program  (usually  the  APL 
directory)  is  your  current  directory. 

Enter:  APLCOM  if  you  are  using  APL  version  5.1  or 
Enter:  A  if  you  are  using  APL  version  6.4. 

These  commands  call  up  batch  files  specifically  for  each  version.     A  batch  file 
contains  a  set  of  commands  that  must  be  executed  before  it  runs  the  APL  program 
in  order  to  have  the  software  matched  with  the  specific  hardware. 
The  system  will  respond  as  follows: 


APL*PLU3  PC  APPLICATION  DEVELOPMENT  SYSTEM 


Version  5.1       Serial   Number  134257         Copyright  1985,    ST3C,    Inc. 
All  rights  reserved.     Unauthorized  reproduction  of  this  software 
is  prohibited  and  violates  U.S.    Copyright  Lays.     APL*PLU3  is  a 
registered  trademark  and  service  marie  of  STSC,    Inc. 

CLEAR  US     11/20/1387  14:31:54 

3.  Place  the  desired  DSP  library  disk  in  the  default  disk  drive,  as  for  example, 
the  disk  containing  the  UTILITY  workspace. 

Enter:  )LOAD  UTILITY 

If  you  are  using  a  hard  disk  and  want  to  load  the  workspace  from  drive  A, 
enter:  )LOAD  0  UTILITY.  The  "0"  (zero)  stands  for  drive  A,  "1"  stands  for 
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drive  B,  and  "2"  stands  for  drive  C  (hard  disk).  The  system  will  respond  as 
follows: 


) LOAD  0  UTILITY 
0  UTILITY  SAVED  12/13/1987  15:12:28 


NAVAL  POSTGRADUATE  SCHOOL  MONTEREY .  CALIFORNIA 

UTILITIES  WORKSPACE 


THESIS  BY 
Y.KATZIR.I.A.F. 
ADVISOR: PROF.  C.W.  THERRIEN. 
VERSION  1.0  SEPTEMBER  1987 


At  this  point  you  can  use  the  UTILITY  workspace  as  described  in  Appendices 
B,  C,  and  D. 

4.  Before  terminating  the  session  save  your  results  by  entering: 

)SAVE  UTILITY  or  )SAVE  0  WSID  (WSID  means  your  file  name) 

the  system  will  respond: 

)SAVE 
0  UTILITY  SAVED    11/19/1987  00:      24:06 

5.  To  terminate  APL  enter:  )OFF  and  you  will  be  automatically  returned  to 
DOS. 

D.    USING  STATGRAPHICS  WITH  APL 

Follow  the  instructions  below  to  start  the  session: 

1.  Start  the  computer  under  DOS.  If  a  graphics  printer  is  attached  be  sure  that 
the  printer  is  set  up  for  graphics  printing.  (DOS  has  to  run  graphics.com  when 
system  boots.) 

2.  Load  APL  as  described  in  the  previous  section  (C.b.) 

3.  After  the  APL*PLUS  PC  system  is  loaded,  replace  the  APL  disk  with  the 
STATGRAPHICS  start-up  disk  and  enter: 

)LOAD  STATGRAF 
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WARNING:  When  using  the  STATGRAF  workspace  use  only  the  )PCOPY 
command  to  enter  functions  or  variables  into  the  workspace. 


During  the  loading  procedure,  the  system  will  display  a  copyright  banner  and 
the  following  message: 
"SYSTEM  INITIALIZING.    PLEASE  BE  PATIENT.    THIS  WILL  TAKE  A  FEW  MOMENTS." 

When  initialization  is  complete,  the  system  will  prompt  for  questions  that  are 
related  to  the  computer  configuration.  Select  the  appropriate  answers  (See  Ref. 
3,  pg.  3-1  to  3-3  for  more  details.)  The  system  will  ask  if  you  want  to  save  the 
settings  previously  entered  in  order  to  have  an  automatic  logon  the  next  time  you 
start  STATGRAPHICS.  The  Main  Menu  will  then  be  displayed  as  in  Figure  A.l. 


! STATGRAPHICS  Statistical  Graphics  Systemi 


DATA  MANAGEMENT  AMD  SYSTEM  UTILITIES 
A.  Data  Management 
3.  System  Environment 
C.  Report  Writer  and  Graphics  Replay 

0.  Plotter  Interface 

PLOTTING  AND  DESCRIPTIVE  STATISTICS 

E.  Plotting  Functions 

F.  Descriptive  Methods 

G.  Estimation  and  Testing 
H.  Distribution  Functions 

1.  Exploratory  Data  Analysis 

ANOVA  AND  REGRESSION  ANALYSIS 
J.  Analysis  of  Variance 
K.  Regression  Analysis 


TIME  SERIES  PROCEDURES 
L.  Forecasting 
M.  Quality  Control 
N.  Smoothing 
0.  Time  Series  Analysis 

ADVANCED  PROCEDURES 

P.  Categorical  Data  Analysis 

Q.  Multivariate  Methods 

R.  Nonparametric  Methods 

S.  Sampling 

T.  Experimental  Design 

MATHEMATICAL  AND  USER  PROCEDURES 
U.  Mathematical  Functions 
V.  Supplementary  Operations 


Use  cursor  keys  to  highlight  desired  section.   Then  press  ENTER. 
lHelp   2Edit   3Savscr  4Prtscr  5       5Go     TVars   8Cmd    9Review  lOQuit 
INPUT     11/18/87   13:08  STATGRAPHICS  Vers.  2.1 


Figure  A.l.    Statgraphics  Main  Menu. 
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4.  When  you  are  in  the  Main  Menu,  press  the  ESC  key  to  enter  APL  immediate 
execution  mode.  Execute  the  function.  Enter:  RESTART  to  return  to  the 
menu-control  mode. 

Since  the  complete  set  of  functions  in  the  STATGRAPHICS  system  exceeds 
the  memory  space  available  on  most  personal  computers,  you  should  load  only 
those  functions  and  variables  that  are  necessary  to  perform  the  desired  tasks.  For 
DSP  applications  the  major  procedures  used  are  described  in  the  following  section. 

D.    STATGRAPHICS  MAJOR  PROCEDURES  USED 

After  starting  STATGRAPHICS,  the  first  display  is  the  Main  Menu  (See  Fig- 
ure A.l.) 

The  Main  Menu  contains  22  sections  under  six  major  categories. 

The  Main  Menu  gives  you  access  to  all  the  data  analysis  and  graphics  proce- 
dures. 

A  highlighted  cursor  bar  is  displayed  on  the  Main  Menu  screen.  In  order  to 
select  any  desired  Menu  Procedure,  move  the  cursor  using  the  arrow  keys  and  press: 
ENTER.  STATGRAPHICS  procedures  are  stored  in  four  disks.  When  using  floppy 
disks  carefully  follow  the  instruction  on  the  screen. 

The  following  is  a  short  description  of  the  most  frequently  used  procedures  in 
this  DSP  software  package: 

1.    System  Environments  [3] 

The  STATGRAPHICS  system  environment  section  is  under  category: 

"DATA  MANAGEMENT  AND  SYSTEM  UTILITIES" 
This  section  allows  you  to  modify  the  STATGRAPHICS  system  environ- 
ment to  fit  your  needs.  This  section  has  nine  procedures.  Those  most  often  used 
are: 
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Screen  Options 


Normal  Text 

Foreground  Color  CO-7): 
Background  Color  CO-7): 
Border  Color  CO-IS): 
Intense  CO-1): 
Blinking  CO-1): 

Active  Text 


Foreground  Color  CO-7) 
Background  Color  CO-7) 
Intense  CO-1): 
Blinking  CO-1): 


Graphics  Mode 

7  *  Resolution  Cl=Low.2=High) :   2 

1  Background  Color  CO-15):     0 

1  Palette  Number  CO-1):       1 

0 

0  *  not  supported  on  all  boards 

Highlighted  Text 


0  Foreground  Color  CO-7):  0 

5  Background  Color  CO-7):  2 

0  Intense  CO-1):  0 

0  Blinking  CO-1):  0 


Figure  A. 2.    The  Screen  Options  Menu. 

1.  Screen  option:  Allows  you  to  change  colors,  text,  and  properties  of  the  graphics 
screen.  (See  Figure  A.2). 

2.  Graphic  options:  Allows  you  to  change  sizes,  orientation,  and  positions  of  the 
graphics  screen.  (See  Figure  A. 3.) 

Graphics  Options 

Line  colors:  123456739ABCDEF1234S 

Point  colors:  56789ABCDEF123451888 

Paint  colors:  39ABCDEF123451234567 

Text  color:  1 


Axis  Color 

Ticmarks 

i  Text 

1 

Size 

Frame:   Yes 

Color      Viewpoint 

X 

1 

In 

1 

X      600 

Y 

1 

In 

1 

1 

Grid:   Both 

1      Y    -1600 

Z 

1 

In 

1 

1 

Border:   None 

1      Z      300 

Titlel 

1 

1 

Title2 

1 

1 

Graphics  Cursor: 
Pointsize:  .005 

:  Yes  Graphics:  Yes 
Ticmark  length:  .040 

Window 

:al 

Origin 

Width 

Split:  1 

Printing  Zoom 

Position:  1 

Horizoni 

.200 

.750 

Factors 

VerticaJ 

.200 

.640 

Vertical:  1 

Horizontal:  1 

Page  Orientation:  Vertical 

Figure  A.3.    Graphics  Options  Menu. 

2.    Plotting  Functions  [3] 

Plotting  functions  are  found  on  the  Main  Menu  under  the  category: 
"PLOTTING  AND  DESCRIPTIVE  STATISTICS" 
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STATGRAPHICS  provides  seven  plotting  procedures  you  can  use  to  pro- 
duce a  variety  of  plots.  The  following  is  an  example  of  using  the  X-Y  Time  Plot 
procedure. 

1.  Start  in  the  APL  environment  and  define  the  variables  X  and  Y  as  follows: 

0  STATGRAF  APL  Ins 
You  are  now  under  control  of  the  APL  interpreter. 
To  resume  menu  control,  type: 

RESTART 
*    X«-J.1Q 

X 
123456739    10 
►         Y«-l    2    2.5    3.3    4    5    6    7.5    3.3    9.5 

Y 

1  2    2.5    3.3    4    5    6    7.5    3.3    9.5 

Enter:  RESTART 
Select  the  X-Y  Line  and  Scatterplots  procedures  from  the  Main  Menu. 

2.  Fill  in  the  Menu  entry  as  demonstrated  in  Figure  A. 4. 

(USE  CAPITAL  LETTERS) 

X-Y  Line  and  Scatterplots 


X:    X 
Y:    Y 


Point  labels: 

Point  codes: 
Point  colors: 
Std.  err.  X: 
Std.  err.  Y: 


Points:  Yes 
Lines :   Yes 


Log 

Mo 

Mo 


Figure  A.4.    X-Y  Plotting  Menu 

3.  Press  the  F6  key  in  order  to  see  the  plot.  (See  Figure  A. 5.) 

4.  Press  the  F4  key  to  print  the  screen,  (if  graphics  printer  is  attached) 

5.  Press  the  F5  key  to  make  any  changes  in  the  plot  (including  titles). 

6.  Press  the  F10  key  to  leave  the  screen  and  to  return  to  the  menu. 
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Figure  A.5.    The  Plot  of  Y  vs  X. 

3.    Response  Surface  Plotting  [3] 

This  procedure  produces  3-D  surface  plots  and  contour  plots,  for  previ- 
ously generated  two  dimensional  functions.  The  procedure  is  under  section:  "Ex- 
perimental Design"  on  the  main  menu  under 

"ADVANCED  PROCEDURES" 
The  following  is  an  example  of  how  to  use  this  procedure: 

1.  Start  with  a  definition  of  a  matrix  of  values  for  a  2-D  function  in  the  APL 
environment  as  in  the  following  example: 
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2.  Restart  STATGRAF  and  select  the  section  "experimental  design."    The  fol- 
lowing will  be  displayed. 


EXPERIMENTAL  DESIGN 


1.  Full  and  Fractional  Factorials 

2.  Central  Composite  Designs 

3.  Alias  Structure 

4.  Response  Surface  Plotting 


3.  Select  "Response  Surface  Plotting"  and  make  the  entry  in  the  menu  as  de- 
scribed in  Figure  A.6. 


Response  Surface  Plotting 


-»  X-axis  minimum:  0 

■*  X-axis  maximum:  16 

■*  Number  of  x-axis  intervals:  15 

■*  Function  type:  Matrix  input 


-»  Y-axis  minimum:  0 

*  Y-axis  maximum:  16 

*  Number  of  y-axis  intervals:  15 


Figure  A. 6.    Surface  Plotting  Menu. 


4.  Press  the  F6  key  and  fill  in  the  entries  as  shown  in  the  example  of  Figure  A. 7. 


Response  Surface   Plotting 


X-axis  minimum:    0 

X-axis  maximum:    16 

Number  of  x-axis   intervals:    15 

Function  type:    Matrix  input 


Y-axis   minimum:    0 

Y-axis  maximum:    16 

Number  of  y-axis   intervals:    15 


Grid  matrix:    XX 


Figure  A. 7.    Selecting  the  Function  Definition. 
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Surface  Plotting  Options 


THE  MATRIX  XX  CEXAMPLE) 


Top  title: 
C2  lines) 
X-Axis  Title:  X 
Y-Axis  Title:  Y 

Z-Axis  Title:  Z 


Plot  type:  Surface 

For  surface  plots: 

X-axis  resolution:  100 
Y-axis  resolution:  100 
Lines  parallel  to:  Both  axes 


For  contour  plot: 
Type:  Lines 
Contour  levels: 

or 
Number  of  Divisions: 


10  Minimum:  1 


Figure  A.8.    Surface  Plot  Entry. 


X-axis  skip  factor:  1 
Y-axis  skip  facror:  1 
Hidden  line  removal:  None 


Maximum:  256 


THEHftTSlXXX(£vAMPT^ 


16 


Figure  A.9.    Surface  2-D  Plot. 

5.  Press  the  F6  key  and  enter  the  values  for  the  surface  or  contour  plot  as  shown 
in  the  example  of  Figure  A.8. 

6.  Press  the  F6  key  to  execute  the  plot,  the  F4  key  to  print  the  plot,  and  the  F5 
key  to  edit  the  plot  (See  Figue  A. 9). 

7.  Press  the  F10  key  to  leave  this  screen  and  return  to  the  menu. 
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WARNING:  Since  the  number  of  points  in  3-D  plots  may  exceed  the  memory 
space  available,  do  not  select  hidden  or  partially  hidden  lines  in  complicated 
plots.  Select  the  entry  "points"  instead  of  entry  "lines"  when  drawing  compli- 
cated contour  plots. 
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APPENDIX  B 
DIGITAL  SIGNAL  PROCESSING  LIBRARY 

Listed  below  are  functions'  names  and  brief  discussions  of  the  various  computer 
algorithms  used  to  accomplish  this  thesis  research.  All  programs  were  written  in 
APL  by  the  author  unless  otherwise  noted.  Program  listings  and  samples  of  runs 
are  given  in  Appendix  C.  For  each  function  its  name  and  the  page  containing  its 
listing  in  Appendix  C  are  given.  Then  the  calling  sequence  is  given,  followed  by  a 
brief  description. 

A.    UTILITY  WORKSPACE 

1.  NORD  (Pg.  68) 

Z<-N  NORD   S 

This  function  generates  a  normally  distributed  random  vector  using  the 
function  UNRD. 

The  function  inputs  are: 

a.  N:  Number  of  elements  in  the  vector 

b.  S:  Vector  with  two  arguments: 
*[Sl]:  The  required  mean 

*[S2]:  The  required  standard  deviation 

The  function's  source  is  the  software  package  GRAFSTAT. 

2.  UNRD  (Pg.  68) 

R<-  N  UNRD   B 
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This  function  generates  a  uniformly  distributed  random  vector.  The  inputs 
are: 

a.  N:  Numbers  of  elements  in  the  vector 

b.  B:  Required  mean 

3.    PUTDATA  CPg.  70) 

F  PUTDATA    'A:\WORK\NAME. EXT' 

This  function  transfers  data  (vector,  array,  etc.)  in  free  format  from  an 
APL  workspace  into  an  IBM-PC  DOS  file.  The  function  uses  APL  NATIVE  files 
[7].  A  bell  sounds  when  the  data  transfer  has  been  completed.  The  function  inputs 
are: 

a.  F:  The  data  to  be  transferred 

b.  'A:\WORK\NAME.EXT':  The  designated  native  file  or  path:  "A : "  Stands  for 
the  drive  type,  "WORK"  stands  for  the  directory,  and  "EXT"  stands  for  extension 
of  the  DOS  file. 


WARNING:  The  user  is  not  permitted  to  designate  a  file  that  already  exists 
as  a  DOS  NATIVE  file. 


4.    GETDATA  (Pg.  72) 

F<-GETDATA    'A:\WORK\NAME. EXT' 

This  function  transfers  data  (vector,  array,  etc.)  from  an  IBM-PC  DOS 
file  to  an  APL  workspace.  The  function  uses  APL  NATIVE  files  [7].  A  bell  sounds 
when  the  data  transfer  has  been  completed. 

1 A :  \W0RK\NAME .  EXT ' :  The  native  file's  name  or  path  has  to  be  transferred: 
"A : "  stands  for  the  drive  type,  "WORK"  stands  for  the  directory,  and  "EXT"  stands 
for  extension  of  the  DOS  file. 


47 


B.    EC3400  WORKSPACE  [4] 
1.    DIGFREQ  (Eg.  75) 

TC«-FS  DIGFREq  FC 

This  function  calculates  the  digital  frequency  using  the  equation: 

e  =  uT  =  2izf-l/fa  (B.l) 

f  =  The  analog  frequency 
fs  =  The  sampling  frequency 
The  inputs  are: 

a.  FS:  The  sampling  frequency  [Hz] 

b.  FC:  The  analog  frequency  [Hz] 

The  output  is: 
TC:  The  digital  frequency  [rad.] 

2.    LPCOEFF  (Pg.  76) 

hLP<-n  LPCOEFF  TC 

This  function  calculates  the  21  +  1  !owpass  FIR  filter  coefficients  using  the 
Frequency  Transformation  and  the  equation: 

hLP(n)  =  V^rfn  -  I)  sin([n  -  J]0C),  n  =  0, 1, 2, . . . ,  21  (B.2) 

where  k  =  1  and  0C  =  The  digital  cutoff  frequency 

The  calculated  coefficients  are  for  a  causal  nonrecursive  filter.  In  order  to 
generate  noncausal  coefficients  use  the  function  SHIFT.  The  inputs  are: 

a.  n:  Number  of  required  coefficients  (Use  an  odd  number.)  By  choosing  a  final 
number  of  coefficients  you  automatically  define  a  rectangular  window. 

b.  TC:  The  digital  cutoff  frequency  [rad.] 
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3.  HPCOEFF  (Pg.  77) 

hHP<-n  HPCOEFF  TC 

This  function  calculates  the  21  +  1  highpass  FIR  filter  coefficients  using 
the  Frequency  Transformation  method  and  the  equation: 

hHp(n)  =  (-l)(n-7)/*LP(n),  n  =  0, 1,2, . . . ,  21  (5.3) 

The  calculated  coefficients  are  for  a  causal  nonrecursive  (linear  phase) 
filter.  In  order  to  generate  noncausal  coefficients  (zero  phase  filter)  use  the  function 
SHIFT.  The  coefficients  are  normalized.  The  inputs  are: 

a.  n:  The  number  of  the  required  coefficients  (Use  an  odd  number.)  By  choosing 
a  final  number  of  coefficients  you  automatically  define  a  rectangular  window. 

b.  TC:  The  digital  cutoff  frequency  [rad.] 

4.  BPCOEFF  (Pg.  78) 

hBP<-n  BPCOEFF  TUL 

This  function  calculates  the  21  +  1  bandpass  FIR  filter  coefficients  using 
the  Frequency  Transformation  method  and  the  equations: 

hBP{n)  =  [2  cos  [(n  -  I)0O]]  hLP(n),  n  =  0, 1, 2, . . .  ,2/  (5.4) 

0c  =  (0u-9e)/2 
0o  =  {9u  +  6t)/2 

9U  =  The  upper  digital  frequency 

$t  —  The  lower  digital  frequency 

Bc  =  The  digital  cutoff  frequency 

The  calculated  coefficients  are  for  a  causal  nonrecursive  filter.  In  order 
to  generate  noncausal  coefficients  use  the  function  SHIFT.  The  coefficients  are 
normalized.  The  inputs  to  the  function  are: 
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a.  n:  The  number  of  the  required  coefficients  (Use  an  odd  number.)  By  choosing 
a  final  number  of  coefficients  you  automatically  define  a  rectangular  window. 

b.  TUL:  Vector  of  two  arguments 

TUL  [1]:  The  upper  frequency  [rad.] 
TUL  [2]:  The  lower  frequency  [rad.] 

5.    BSCOEFF  (Pg.  79) 

hBS*-n  BSCOEFF  TUL 

This  function  calculates  the  21  +  1  bandstop  FIR  filter  coefficients  using 
Frequency  Transformation  and  the  equation: 

hBS(0)  =  1  -  hBp(Q) 

(B.5) 
hBs(n)  =  -hBp(n),n  =  1,2, ...  ,27 

Bc  =  (0U  -  0e)/2 

0o  =  (On  +  9e)/2 

9U  =  The  upper  digital  frequency 

9f  =  The  lower  digital  frequency 

9C  =  The  digital  cutoff  frequency 

The  calculated  nonrecursive  coefficients  are  for  a  causal  nonrecursive  filter. 
In  order  to  generate  noncausal  coefficients  use  the  function  SHIFT.  The  coefficients 
are  normalized.  The  inputs  to  the  function  are: 

a.  n:  The  number  of  the  required  coefficients  (Use  an  odd  number.)  By  choosing 
a  final  number  of  coefficients  you  automatically  define  a  rectangular  window. 

b.  TUL:  Vector  of  two  arguments 

TUL  [1]:  The  upper  frequency  [rad.] 
TUL  [2]:  The  lower  frequency  [rad.] 
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6.    Hamming  (Pg.  80) 

W<-HAMMING  N 

This  function  generates  an  27  +  1  point  Hamming  window  using  the  equa- 
tion: 


W(n)  =  0.54  +  0.46  cos 


7r(n  —  7) 


,n  =  0,l,2,...,2J  (B.6) 


I 

The  window  starts  at  n  =  0  and  has  27  -f- 1  samples.  The  input  to  the  function  is 
N:  The  number  of  samples. 

7.    HANNING  (Pg.  81) 

Wf-HANNING  N 

This  function  generates  an  27  +  1  point  Hanning  (von  Hann)  window  using 
the  equation: 


W(n)  =  0.5  +  0.5cos 


7r(n  —  7) 


,n  =  0,l,2,...,2I  {B.l) 


The  window  starts  at  n  =  0  and  has  27+1  samples.    The  input  to  the 
function  is: 
N:  The  number  of  samples. 

8.    BARTLETT  (Pg.  82) 

W<-BARTLETT  N 

This  function  generates  an  27+1  point  Bartlett  (triangular)  window  using 

the  equation: 

f  1  -  i(n  -  7)     ,0  <  n  <  27 
W(n)  =  I  (B.S) 

{  0  ,  otherwise 

The  window  starts  at  n  =  0  and  has  27+1  samples.    The  input  to  the 
function  is: 
N:  The  number  of  samples. 
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9.    FREQRES  (Pg.  83) 

HW<-Z  FREQRES   h 

This  function  generates  the  magnitude  and  phase  of  the  frequency  response 
from  the  unit  impulse  response.  The  functions:  FFT ,  XMAGN ,  XPHASE  are  used.  In 
order  to  achieve  good  resolution  the  coefficients  are  zero  padded  to  'Z'  samples. 


WARNING:  The  'Z'  needs  to  be  a  power  of  2. 


The  inputs  to  the  function  are: 

a.  Z:  The  total  number  of  samples  including  the  zero  padding 

b.  h:  The  impulse  response  coefficients 
The  outputs  are: 

a.  HW:  The  frequency  response  magnitude  (Z/2  samples) 

b.  PW:  The  frequency  response  response  phase  (Z/2  samples)  [rad.] 

10.    IDEALF  (Pg.  84) 

H<-TYPE  IDEALF  TC 

This  function  produces  the  frequency  response  of  various  ideal  filters.  The 
inputs  are: 

a.  'TYPE':  The  tvpe  of  the  filter  to  be  generated: 
'LP':  Lowpass 
'HP':  Highpass 
'BP':  Bandpass 
'BS':  Bandstop 
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b.  TC:  A  vector  containing  the  following: 

1.  In  the  case  of  LP  or  HP: 

TC  [1]  =  The  number  of  samples 

TC  [2]  =  The  cutoff  digital  frequency 

2.  In  a  case  of  BP  or  BS: 

TC  [1]  =  The  number  of  samples 

TC  [2]  =  The  upper  digital  frequency 

TC  [3]  =  The  lower  digital  frequency 

The  frequencies  should  be  normalized  to  lie  in  the  interval  0  to  it.  [rad.] 
11.    FCOEFF  (Pg.  86) 

h<-n  FCOEFF  HW 

This  function  calculates  the  27+1  point  filter  coefficients  using  IFFT  (The 
sampling  method.) 

The  coefficients  are  for  a  causal  nonrecursive  filter.    The  inputs  to  the 
function  are: 

a.  n:  The  number  of  the  required  coefficients  (Use  an  odd  number.)  By  choosing 
a  final  number  of  coefficients  you  automatically  define  a  rectangular  window. 

b.  HW:  The  ideal  frequency  response. 
C.    EC3410  WORKSPACE  [5,  8] 

1.    SACF  (Pg.  88) 

R<-LAG  SACF  DAT 

This  function  computes  the  sample  autocorrelation  function  (unbiased) 
using  the  following  equation: 

AT-1-|*I 

^  =  1/^      Y,     Yn+\k\+Yn,     0<k<N-l  (B.9) 

n=0 
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The  inputs  are: 

a.  DAT:  The  data  block  (one- dimensional  vector) 

b.  LAG:  Number  of  lags  to  be  calculated  (the  limit  is  300  lags) 

2.  MEAN  (Pg.  90) 

M«-MEAN  F 

This  function  calculates  the  sample  mean  of  a  data  set  using  the  following 
equation: 

N 

m  =  l/Nj2Y(i)  (B-10) 

The  input  is: 
F:  The  data  block 

3.  COV  (Pg.  91) 
K<-C0V  R 

This  function  generates  a  Toeplitz  covariance  matrix  using  covariance  func- 
tion. The  input  is: 
R:  The  covariance  function  values  (vector) 

4.  LPFP  (Pg.  92) 
a«-P  LPFP  DATA 

This  function  generates  the  pth  order  Linear  Predictive  Filter  parameters 
and  calculates  the  Prediction  Error  Variance  by  solving  the  Normal  Equation: 

R(0)        R(l)     ...     R{p-l) 


R{\) 

IR(P-I)      . 
The  inputs  are: 
a.  P:  The  coefficients  order 


R(0) 


1   1 

\a<] 

a\ 

= 

0 

-ap-i . 

.  0  . 

(5.11) 
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b.  DATA:  The  data  block 
The  outputs  are: 

a.  a:  The  filter  coefficients 

b.  PRER:  The  estimate  prediction  error  variance 
The  following  functions  are  used: 

a.  SACF:  Calculates  the  sample  autocorrelation  of  the  data  block 

b.  COV:  Generates  the  Toeplitz  matrix  of  the  covariance  function 

5.    CC  (Pg.  93) 
Y<-x  CC  h 

This  function  performs  one-dimensional  circular  convolution  between  two 
inputs  (h(n)  0i(n)).  x  and  h  are  the  arguments  to  be  convolved. 


WARNING:  x  and  h  need  to  be  the  same  size  and  one- dimensional, 
otherwise  an  error  message  will  be  generated  and  the  convolution  will  not 
start. 


6.    FCC  (Tg.  94) 

Y<-x  FCC  h 

This  function  does  one- dimensional  fast  circular  convolution  between  two 
arguments.  The  fast  operation  is  established  by  a  shifted  matrix  multiplication. 
x  and  h  are  the  arguments  to  be  convolved.  Due  to  memory  size  this  function 
is  limited  to  arguments  of  less  than  80  numbers.  If  the  argument  is  too  long  a 
printout  message  advises  the  user  to  use  regular  circular  convolution  function. 


WARNING:  x  and  h  need  to  be  the  same  size  and  one-dimensional, 
otherwise  an  error  message  will  be  generated  and  the  convolution  will  not 
start. 
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7.    LCV  (Pg.  95) 

Y«-xl  LCV  hi 

This  function  does  one  dimensional  linear  convolution  between  two  argu- 
ments by  using  the  following  equation: 


oo 


xi(n)  *  hi(n)  =     y_]    xi(m)hi(n—m)  (5.12) 


m=  —  oo 


xl  and  h\  are  the  arguments  to  be  convolved. 

The  function  uses  the  CC  function  (Circular  Convolution)  and  the  linear 
convolution  is  produced  by  zero  padding  the  arguments. 

8.  FLCV  (Pg.  96) 
Y<-xl  FLCV  hi 

This  function  does  one- dimensional  linear  convolution  (same  as  LCV  func- 
tion) by  using  the  function  FCC  in  order  to  achieve  a  faster  computation  time. 

9.  FFT  (Pg.  97) 
Z*-FFT  X 

This  function  calculates  the  one-dimensional  Discrete  Fourier  Transform 
according  to  the  definition 

N-l 

X(k)  =  Y,  x(n)e-*2ltlN)nk,k  =  0, 1, . . . ,N -  1  (5.13) 

n=0 

The  input  to  the  function: 

X:  Can  be  either  real  vector  of  length  N  or  a  two  by  N  matrix  (the  1st  row  is  the 
real  part  and  the  2nd  row  is  the  imaginary  part).  N  must  be  a  power  of  two. 

The  output  is: 

Z:  A  complex  two  by  N  matrix  (real  and  imaginary) 

The  function  was  written  by  Professor  Paul  Penfield,  Jr.  of  M.I.T. 
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10.    IFFT  (Pg.  98) 
Z<-IFFT  X 

This  function  calculates  the  one-dimensional  Inverse  Discrete  Fourier  Trans- 
form according  to  the  definition: 

"JV-l 


x(n)  =  l/N 


£  *(*> 


-j(2n/N)nk 


Jt=0 


,n  =  0,l,...,iV-l 


(5.14) 


The  input  to  the  function  is: 


X:  Can  be  either  a  real  vector  of  length  N  or  a  two  by  N  matrix  (the  1st  row  is 
the  real  part  and  the  2nd  row  is  the  imaginary  part).  N  must  be  a  power  of 
two. 

The  output  is: 

Z:  A  complex  two  by  N  matrix  (real  and  imaginary) 

The  function  was  written  by  Professor  Paul  Penfield,  Jr.  of  M.I.T. 
11.    PSE  (Pg.  100) 
INW<-ZP  PSE  DATA 

This  function  computes  the  power  spectrum  estimate  (Periodogram)  using 
the  equation: 


Sx  (e?0)  =  IN (6)  =  l/N 


X  (*>') 


(£.15) 


The  function  uses  the  FFT  APL  function  and  the  inputs  are: 

a.  DATA:  The  data  block  (one  dimensional  vector) 

b.  ZP:  The  length  of  the  zero  padding  needed  for  a  smooth  plot 


WARNING:  The  DATA  plus  ZP  need  to  be  in  length  which  is  a  power  of  two 
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12.    PSEB  (Pg.  100) 

INBW^-SEG   PSEB   DATA 

This  function  computes  the  power  spectrum  estimate  with  rectangular 
nonoverlapping  window  using  the  equation: 

2 


IN'(i\6)  =  l/N 


X&  (e>*) 


M 


(£.16) 


B{e)  =  i/M^iN'^\e) 


t=i 


N:  The  length  of  the  data  set 
M:  The  length  of  a  segment 

Each  segment  has  N'  Points:  N  =  N'  x  M. 

This  function  uses  the  FFT  and  the  PSE  APL  functions  and  the  inputs  are: 

a.  SEG:  The  length  of  the  segments 

b.  DATA:  The  data  block  (one- dimensional  vector) 

The  function  automatically  zero  pads  the  segment  up  to  256  points. 


WARNING:  In  case  of  error  due  to  memory  size,  use  fewer  segments 


13.    XSQRT  (Pg.  102) 

XZ<-XSQRT  XR 

This  function  calculates  the  complex  square  root  of  a  complex  number. 

XZ  =  XR1/2  (B.17) 
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14.  XTIME  (Pg.  102) 
XZ*-XL  XTIME  XR 

This  function  multiplies  two  complex  numbers 

XZ  =XLxXR  (5.18) 

15.  XMAGN  (Pg.  103) 

Z<-XMAGN  XR 

This  function  computes  the  complex  magnitude  of  a  number,  the  result  is 
a  real  number 

Z=\XR\  (5.19) 

16.  XPHAS  (Pg.  103) 

Z<-XPHAS  XR 

This  function  computes  the  phase  of  the  complex  number  XR,  and  returns 
a  real  number. 

Z  =<£  XR  (5.20) 

17.  XCONJ  (Pg.  104) 

XZ<-XC0NJ  XR 

This  function  returns  the  complex  congugate  of  the  input. 

18.  XEPO  (Pg.  104) 
XZ<-XEP0   XR 

This  function  computes  the  complex  exponential 

XZ  =  eXR  (5.21) 


19.    PI  (Pg.  105) 
Y«-PI 
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This  function  returns  the  it  value. 

20.  COS  (Pg.  105) 
Y<-C0S  X 

This  function  returns  the  cosine  of  the  input. 

21.  SIN  (Pg.  105) 
Y<-SIN  X 

This  function  returns  the  sine  of  the  input. 

D.    EC4440  WORKSPACE  [9] 
1.    CC2D  (Pg.  107) 

Y<-X  CC2D  H 

This  function  calculates  two  dimensional  circular  convolution  using  the 
equation: 

Y  =  h   *  *  x 

Nr-l  N2-l 

Y{nl,n2)=    Y^    X]  ^(ni'n2)  x  [((ni  ~mi))ivi((u2  -m2))Ni]  (£.22) 

nl  =  0   n2=0 

Where  subscript  N\  and  N2  denotes  arithmetic  modulo  ivi  and  N2.    x  and  h  are 
the  arguments  to  be  convolved. 


WARNING:  X  and  H  need  to  be  two-dimensional  with  the  same  size,  other- 
wise an  error  message  will  be  generated  and  the  convolution  will  not  start. 


2.    FFT2D  (Pg.  109) 

Y^-FFT2D  X 

This  function  calculates  the  two-dimensional  DFT  using  the  one-  dimen- 
sional FFT  function.  If  the  input  is  real,  the  function  converts  it  to  a  complex  array. 
If  the  output  array  is  real  valued  the  second  plane  is  eliminated  automatically. 
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WARNING:  The  input  needs  to  be  a  length  which  is  a  power  of  two. 


3.  IFFT2D  (Pg.  110) 
Y<-IFFT2D  X 

This  function  calculates  the  two-dimensional  Inverse  DFT  using  one-dimensional 
IFFT  function  (the  same  explanation  as  FFT2D). 

4.  SHIFT  (Pg.  114) 
Y<-SHIFT  X 

This  function  shifts  periodic  data  to  an  interval  starting  at  the  origin. 

5.  UNSHIFT  fPg.  114) 
Y<-UNSHIFT  X 

This  function  shifts  periodic  data  to  an  interval  centered  at  the  origin. 

6.  SHIFT2D  (Pg.  115) 
Y<-SHIFT2D  X 

This  function  shifts  two-dimensional  periodic  data  to  an  interval  starting 
at  the  origin 

7.  UNSHIFT2D  (Pg.  115) 
Y<-UNSHIFT2D  X 

This  function  shifts  two-dimensional  periodic  data  to  an  interval  centered 
at  the  origin.  Each  of  the  following  functions  generates  a  two-dimensional  window 
by  using  a  one-dimensional  window.  Two  types  of  supports  are  used: 

a.  A  rectangular  region  of  support  which  is  achieved  by  computing  the  outer 
product  of  two  1-D  windows: 

WR(nun2)  =  W1(n1)W2(n2)  (£.23) 
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b.  A  circular  region  of  support  which  is  formed  by  sampling  a  circularly  rotated, 
1-D,  continuous  function: 


Wc{nx,n2)  =  w(jn\+ni\  (5.24) 

8.  ILPFILT  (Pq.  117) 
W<-R  ILPFILT  P 

This  function  generates  a  circular  based  rectangular  window  using  the  1-D 
equation: 

This  can  also  be  used  for  defining  the  frequency  response  of  an  ideal  circular  sym- 
metrical lowpass  filter.  The  inputs  are: 

a.  R:  The  window  base  radius 

b.  P:  The  dimension  of  the  square  matrix  in  which  this  window  is  defined. 

9.  HAMMINGR  (Pg.  119) 

WR*-II   HAMMINGR  N 

This  function  generates  a  2-D  Hamming  rectangular  based  window  using 
the  1-D  equation: 

W(n)  =  I  °'54  "  °-46  cos  (W<0     ,  M  <  I  (£.26) 

\  0  ,  otherwise 

The  inputs  are: 

a.  II:  The  window  base  dimension  (one-side) 

b.  N:  The  dimension  of  the  square  matrix  in  which  this  window  is  defined. 


62 


10.  HAMMINGC  (Pg.  121) 

WC*-II   HAMMINGC  N 

This  function  generates  a  2-D  Hamming  window  with  circular  base  using 
equation  (B.26).  The  inputs  are: 

a.  II:  The  window  base  radius 

b.  N:  The  dimension  of  the  square  matrix  in  which  the  window  is  defined. 

11.  HANNINGR  (Pg.  123) 

WR«-II   HANNINGR  N 

This  function  generates  a  2-D  Hanning  window  with  rectangular  base  using 
the  1-D  equation: 

W(n\  _  J  °-5  I1  +  cos  ( W  *)]     i  M  <  I  (B.27) 

\  0  ,  otherwise 

The  inputs  are: 

a.  II:  The  window  base  dimension  (one  side) 

b.  N:  The  dimension  of  the  square  matrix  in  which  the  window  is  defined. 

12.  HANNINGC  (Pg.  125) 

WC«-II   HANNINGC  N 

This  function  generates  a  2-D  Hanning  window  with  circular  base  using 
equation  (B.27).  The  inputs  are: 

a.  II:  The  window  base  radius 

b.  N:  The  dimension  of  the  square  matrix  which  this  window  is  defined. 

13.  BARTLETTR  (Pg.  127) 

WR«-II   BARTLETTR  N 
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This  function  generates  2-D  Bartlett  (triangular)  window  with  rectangular 
base  using  the  1-D  equation: 


W(n)  =  { 


(  I-  jn     ,  0  <n  <  I 


1  +  jn     , -l<n<0  (5.28) 

0  ,  otherwise 


The  inputs  are: 

a.  II:  The  window  base  radius 

b.  N:  The  dimension  of  the  square  matrix  in  which  the  window  is  defined. 

14.  BARTLETTC  (Pg.  129) 

WC<-II   BARTLETTC  N 

This  function  generates  a  2-D  Bartlett  window  with  rectangular  base  using 
equation  (B.28).  The  inputs  are: 

a.  II:  The  window  base  radius 

b.  N:  The  dimension  of  the  square  matrix  in  which  the  window  is  defined. 

15.  RECTANR  (Pg.  131) 

WR<-II   RECTANR  N 

This  function  generates  2-D  rectangular  window  with  rectangular  base 
using  equation  (B.25).  The  inputs  are: 

a.  II:  The  window  base  dimension  (one  side) 

b.  N:  The  dimension  of  the  square  matrix  in  which  the  window  is  defined. 

16.  PROTFILT  (Pg.  133) 

HPR0T<-PR0TFILT  SIZE 

This  function  produces  a  frequency  response  1-D  lowpass  prototype  filter 
in  order  to  design  a  2-D  FIR  filter  using  a  transformation.  The  function  SHIFT 
has  been  used.  The  input  is: 
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SIZE:  The  size  of  the  prototype  filter  (number  of  samples). 

17.    COEFF  (Pg.  135) 
A«-C0EFF  HP 
This  function  finds  the  coefficients  a(n)  in  the  representation: 

N 


H(u>)  =  2_]  a(n)cos(un) 


n=0 

where: 

a(n)  £  | '£&>       •n  =  °  (B.29) 

v   ;       1 2h(n)     ,n  >  0  v         y 

The  function  IFFT  has  been  used.  The  input  is: 
HP:  The  filter's  1-D  prototype. 

18.  TRANSFNC  (Pg.  137) 
Fvv<-TRANSFNC   SIZE 

This  function  produces  a  transformation  matrix  to  design  a  FIR  filter 
using  McCLELLAN  transformation.  In  this  example  the  following  function  has 
been  used: 

F(<jJi,u>2)  =  1/2 (  —  1  +  costal  +  cosu;2  +  cosu;i  cosu^)  (£.30) 

In  order  to  change  F(u>i ,  uj2  )  to  a  different  transformation  function,  replace 
lines  10,  13,  and  14  in  the  program  with  the  new  choices.  (See  example  in  Appendix 
C  (Pg.  137).  The  input  is: 

SIZE:  The  number  of  samples  of  u\  or  u2. 

19.  CHEB  (Pg.  140) 
Y<-N  CHEB  X 
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This  function  calculates  the  Chebyshev  Polynomials  using  the  following 
recursion: 

Tn(x)  =  2xTn.1(x)  -  Tn_2(x)  (£.31) 

The  inputs  are: 

a.  N:  The  order  of  the  desired  polynomial 

b.  X:  The  argument. 


WARNING:  The  combination  of  high  order  of  Chebyshev  polynomial  and  a 
large  array  to  be  evaluated  will  result  in  a  long  period  of  calculating  time. 


20.    MCCLEL  fPg.  141) 

Hwv«-a  MCCLEL  F 

This  function  computes  the  2-D  frequency  response  using  McCLELLAN 

transformation  according  to  this  equation: 

N 
H{ux ,  u2 )  =  J2  <*(n)Tn  [*X"i ,  w2 )]  (B.32) 

n  =  0 

In  order  to  indicate  the  calculating  status,  the  function  prompts  the  present  cal- 
culating order  of  the  polynomial.  The  inputs  are: 

a.  a:  The  impulse  response  coefficients  of  the  prototype  filter 

b.  F:  The  transformation  matrix  (previously  calculated  using  the  function  TRANSFNC) 
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APPENDIX  C 
FUNCTION  LISTINGS  AND  EXAMPLES  OF  USE 


)LOAD  0  UTILITY 
0  UTILITY  SAVED  12/10/1987  16:37:45 


NAVAL  POSTGRADUATE  SCHOOL 


MONTEREY , CALIFORNIA 


VERSION  1.0 


UTILITIES  WORKSPACE 

THESIS  BY 
Y.KATZIR.I.A.F. 
ADVISOR: PROF.  C.W.  THERRIEN. 

SEPTEMBER  1987 


)FNS 
GETDATA  INFORMATION 


NORD 


PUTDATA  UNRD 
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vNORDcnnv 

CO]  Z<-N  NORD  P;S;I;T 

CI]  fl 

C2]  fl  A  NORMAL  DISTRIBUTION  RANDOM  VECTOR  GENERATOR 

C3]  fl 

C4]    P«-2TP,1 

C5]    Z<-CN)p0 

C6]    I«-l 

[73  F10:T<-2  DNRD  0.5  0.5 

C83  T«-C2xT)-l 

C93  S«-(TC13*2]+TC23*2 

C103  -»F10xt.S£1 

C113  ZCI3«-PC13+PC23xCTC13xCC~2x«S3*S)*0.5D 

[123  ■*F10xUI«-I*13*N 

[13]  -+0 

[14]  fl 

[15]  fl    SOURCE:  GRAFSTAT 

[16]  A    N«-NUMBER  OF  ELEMENTS  IN  THE  VECTOR 

[17]  n    P<-MEAN 

[18]  A    S«-SIGMA  (NOT  SIGMA*2) 

[19]  A 


VUNRDCD]? 

[0]  R«-N  UNRD  B 

[1]  A 

[2  3  B   A  UNIFORM  DISTRIBUTION  RANDOM  VECTOR  GENERATOR 

[33  fl 

[4  3  R«-CB[l]-B[2])  +  CCN?10000)x2xB[2]  5*10000 

[5]  -»0 

[6]  B 

[7]  B   SOURCE:  GRAFSTAT 

[8]  B   N<-NUMBER  OF  ELEMENTS  IN  THE  VECTOR 

[9]  B   B«-THE  LIMITS  C2  NUMBERS) 


[10]  fl 
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GAUSD+-100  NORD  0  1 


GADSD 
1.558274527  "2.444846636  "1.09846582  1.122362631  0.580943715 
"0.2714227577  0.4143999987  "0.9733729641  "1.021759908 
0.3179468935  1.515647616  0.7495965254  "0.5068186108  0.8851652886 
"0.2479112202  "0.7262416263  "0.4452313973  "0.6121927891 
"0.2086518602  0.5618385258  "1.064013245  0.3513946699  1.132830543 
0.1510880125  0.7028398161  "0.05231570544  2.018906718 
0.9240866379  "1.814329439  0.0351494428  "1.806312062  1.028100377 
0.3944470661  0.6391243468  0.8738735133  1.75599265  "0.3197284973 
"0.1366870369  0.6159939661  0.9774365509  "1.115427909 
"0.5496175242  0.03996605725  "2.484219992  1.158618634 
"1.026625814  1.153535117  "0.7858857486  0.6349969194  0.8198888762 
"0.1760331562  0.5627148978  "0.127309154  0.5538766995  "1.09726933 
"0.7311467753  1.404690913  "0.6203249336  0.2370824053 
"1.585608223  "0.4013458241  "0.7707101307  "0.2637646461 
0.9762420248  0.9775308888  1.169573915  0.1594895548  0.5001385714 
"1.05500756  "0.4509335123  1.271139261  0.8989835792  0.438888127 
"1.247429466  0.3241026886  0.3904619568  "0.4052863418 
0.2929672888  2.565225342  "0.4576884119  "1.611439324  "2.669649901 
"0.7594268721  "0.6750622611  "1.171739437  2.032501263 
0.9683033719  0.6699827805  0.4199990777  "2.871065022  1.686296301 
0.02722627758  "0.9020480504  "2.053207747  0.08940644014  2.0866291 
0.3651320115  0.8457720835  "0.1842246182  1.03036774 


UNIFD<-100  UNHD  0  100 


UNIFD 
3.22  "5.78  18.74  "10.04  43.86  91.44  "63.84  "19.84  "18.58  94.66  "51.14 
"46.78  "80.44  "23.98  "67.32  69.64  71.06  "82.22  "14.14  "37.6 
"66.84  70.92  "93.12  "30.38  18.1  92.06  50.94  "90.96  82.12  24.6 

29.6  83.32  "1.1  34.96  "92.3  61.72  "81.84  91.02  "77.52  3.34  "18.2 
26.18  "16.3  "93.02  "75.96  "14.68  66.88  44.58  "75.88  "2.82  "10.68 
"87.08  "34.64  81.86  "61.82  "68.26  48.36  35.22  6.46  "76.02  "84.16 

26.7  74.7  "62.86  "31.94  0.42  "42.96  96.16  33.74  "52.68  39.3 
"51.16  "8.2  "77.86  "60.26  "3.18  "13.76  16.38  64.04  "5.46  "24.86 
"31.52  "58.7  "35.88  "16.94  "47.26  73.98  "87.34  "94.42  "62.04 
"80.84  61.58  "89.98  57.78  "99.9  43.1  60.76  "32.04  "1.38  "89.56 
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?PUTDATACD]V 

C03  P  PUTDATA  »AME;DI0;TIE:FA;NREC;C:r8:FC;SELLJI 

C2]  a  FUNCTION  TO  WRITE  DATA  rvi    nr*^ 

C33  J  INTO  AN  IBM  pX^S.  SleT^Se^^^'  )IN  ****   F0RMAT 

t5]  QI0«-I<-1 

161  ONUNTIE  "1 

C7]  TIE*"1 

C8  3  a  CREATING  A  TIED  FILE 

£9  J  NAME  ONCREATE  TIE 

CIO]  q  WRITING  OUT  THE  SHAPE  OF  F 

C113  FA<-".  io  0  *PF 

C12J  a  APPENDING  THE  SHAPE  OF  F  TO  m  pttt- 

13   cqtclf.otcnl.fadqnappfSd  tS 
c14j  a  reshaping  f         " 

CIS]       NREC«-rCx/pF)-r4 
C16]       01+41  Cx/CpF)}-l 
C17]       F«-CNREC,4)pF 

r??J   C,°Jc"'-DTCNL.FB)ONAPPEND  TIE 

III]    ^^S^"™*  «"-»■««  RECORO 
J25]   FCCCFC«'-')/t.pPC3;,., 

C26]   COTCLF.DTCNL.FOQNAPPEND  TIE 

C27]  a  UNTYING  THE  FILE 

C28]   ONUNTIE  "1 

C29]   BELL*-DTCBEL 

1 3  0  ]   BELL 

C31]   -»o 

C32]  A 

C33]  A  THE  FUNCTION  USES  NATIVES  APL  COMMANDS. 

C35]  A  Y.KATZIR  .I.A.F..  6.1.87 

36   A  ADVISOR:  C.W.  THERRIEN. 
I  J7  J  A 
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GAGSD  POTDATA  'GA0S.DAT* 


)OFP 


C:\>TfPE  GAUS.DAT 


100 
1.558274527 
0.580943715 
-1.021759908 
-0.5068136108 
-0.4452313973 
-1.064013245 

0.7028398161 
-1.314329439  ' 
0.3944470661 
-0.3197234973 
-1.115427909 
1.158613634 
0.6349969194 
-0.127309154 
1.404690913 
-0.4013453241 
0.9775308888 
-1.05500756 

0.438888127 
-0.4052863418 
-1.611439324 
-1.171739437 

0.4199990777 
-0.9020480504 
0.3651320115  0, 


'-2.444846636 

-0.2714227577 
0.3179468935 
0.3351652336 
-0.6121927391 
0.3513946699 
-0.05231570544 
0.0351494423 
0.6391243463 
-0.1366870369 
-0.5496175242 
-1.026625314 
0.3198888762 
0.5538766995 
-0.6203249336 
-0.7707101307 

1.169573915 
-0.4509335123 
-1.247429466 

0.2929672383 

-2.669649901 

2.032501263 

-2.371065022 

-2.053207747 


-1.09846582 
0.4143999987 
1.515647616 
-0.2479112202 
-0.2086513602 
1.132330543 
2.013906713 
-1.306312062 
0.3738735133 
0.6159939661 
0.03996605725 
1.153535117 
-0.1760331562 
-1.09726933 

0.2370824053 
-0.2637646461 
0.1594895548 
1.271139261 
0.3241026886 
2.565225342 
-0.7594268721 
0.9683033719 
1.686296301 
0.08940644014 


3457720835  -0.1342246132  1.0303 


1.122352631 
-0.9733729641 
0.7495965254 
-0.7262416263 
0.5613385258 
0.1510880125 
0.9240866379 
1.023100377 
1.75599265 
0.9774365509 
-2.484219992 
-0.7858857486 
0.5627148973 
-0.7311467753 
-1.535608223 
0.9762420243 
0.5001335714 
0.3989835792 
0.3904619568 
-0.4576384119 
-0.5750622611 
0.6699827305 
0.02722627758 
2.0866291 
6774 
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VGETDATACCHV 

CO]  F«-GETDATA  NAME;0IO;DATA:M:SHAP;I :BELL 

CI]  ft 

C2]  ft  FUNCTION  TO  READ  FREE  FORMATED  DATA  CF]  C VECTOR. ARRAY ,ECT. ) 

C3]  ft  FROM  IBM  P.C.  D.O.S.  FILE  C'NAME'] 

C4]  ft 

C5]    QIO«-I«-0 
C6]    ONUNTIE  "1 

C7]  ft  TYING  THE  FILE 

C8]  NAME  DNTIE  "1 

C9]  ft  FINDING  THE  SIZE  OF  THE  FILE. 

CIO]  M^dNSIZE  ~1 

Cll]  ft  READING  THE  DATA  INTO  THE  WORK  SPACE 

C121  FH3NREAD  -1.82.M 

C13]  FCCF=DAVC45])/T.pF]«-0.AVC253] 

C14]  ft  EXECUTING  THE  DATA  AS  AN  APL  EXPRESSION 

CIS]  F<-*(~CF=0AVC10])vCF=0AVC13])vCF=DAVC26]))/F 

C16]  DNUNTIE  "1 

C17]  LOOP:I<-I  +  l 

C18]  -»(I£10) /ERROR 

C19]  DATA«-HF 

C20]  SHAP«-ITF 

C21]  -»(Cx/SHAP)9«CpDATA))/LOOP 

C22]  F<-SHAPpDATA 

C23]  BELL«-OTCBEL 

C24]  BELL 

C25]  -»0 

C26]  ERROR: 'FILE  INCORECTLY  FORMATED' 

C27]  -»0 

[23]  ft 

C29]  ft  DATA  HAVE  BEEN  CREATED  BY  APL. PL/ I .FORTRAN, EDITOR  ETC. 

C30]  ft 

C31]  ft  NUMBERS  NEED  NOT  HAVE  DECIMAL  POINTS  AND  ARE  DELIMITED  BY  SPACE. 

C32]  ft  FIRST  RECORD  IS  SHAPE  OF  DATA. 

Cj3]  ft 

C34]  ft  Y.KATZIR.I.A.F. .  6.1.87. 

C35]  ft  ADVISOR:  C.W.  THERRIEN. 

C36]  ft 
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GAUSD«-GETDATA  '  GAUS .  DAT ' 


GAUSD 
1.558274527  "2.444846636  "1.09846582  1.122362631  0.580943715  "0.27142275 
0.4143999987  "0.9733729641  "1.021759908  0.3179468935  1.515647616 
0.7495965254  "0.5068186108  0.8851652886  "0.2479112202  "0.726241626 
"0.4452313973  "0.6121927891  "0.2086518602  0.5618385258  "1.06401324 
0.3513946699  1.132830543  0.1510880125  0.7028398161  "0.05231570544 
2.018906718  0.9240866379  "1.814329439  0.0351494428  "1.806312062 
1.028100377  0.3944470661  0.6391243468  0.8738735133  1.75599265 
"0.3197284973  "0.1366870369  0.6159939661  0.9774365509  "1.115427909 
"0.5496175242  0.03996605725  "2.484219992  1.158618634  "1.026625814 
1.153535117  "0.7858857486  0.6349969194  0.8198888762  "0.1760331562 
0.5627148978  "0.127309154  0.5538766995  "1.09726933  "0.7311467753 
1.404690913  "0.6203249336  0.2370824053  "1.585608223  "0.4013458241 
"0.7707101307  "0.2637646461  0.9762420248  0.9775308888  1.169573915 
0.1594895548  0.5001385714  "1.05500756  "0.4509335123  1.271139261 
0.8989835792  0.438888127  "1.247429466  0.3241026886  0.3904619568 
"0.4052863418  0.2929672888  2.565225342  "0.4576884119  "1.611439324 
"2.669649901  "0.7594268721  "0.6750622611  "1.171739437  2.032501263 
0.9683033719  0.6699827805  0.4199990777  "2.871065022  1.686296301 
0.02722627758  "0.9020480504  "2.053207747  0.08940644014  2.0866291 
0.3651320115  0.8457720835  "0.1842246182  1.03036774 
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JLOAD  0  EC3400 
0  EC3400  SAVED  12/12/1987  21:23:25 


NAVAL  POSTGRADUATE  SCHOOL 


MONTEREY . CALIFORNIA 


ECS  -3400  WORKSPACE 


THESIS  BY 
Y.KATZIR.I-A.F. 
ADVISOR: PROF.  C.W.  THERRIEN. 


VERSION  1.0 


NOVEMBER  1987 


)FNS 

BARTLETT         BPCOEFF  BSCOEFF  CC       COS      DIGFREQ  FCC      FCOEFF   FFT 
FLCV     FREQRES  HAMMING  HANNING  HPCOEFF  IDEALF   IFFT     INFORMATION      LCV 
LPCOEFF  MEAN     PI       SHIFT    SIN      UNSHIFT  XCONJ    XEXPO    XMAGN    XPHftS 
XSQRT    XTIME 
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VDIGFREQCCn? 

CO]  TOFS  DIGFREQ  FC 

CI]  fl 

C2]  A   CALCULATING  THE  DIGITAL  FREQUENCY 

C3]  A 

C4]    TC«-C2xPIxFC)+FS 
C5]    -»0 

C6]  A 

C7]  A   FS«-SAMPLING  FREQUENCY  CHz] 

C8]  A   FC*-THE  CORRESPONDING  ANALOG  FREQUNCY  CHz] 

C9]  A   TC:THE  DIGITAL  FREQUENCY  CRAD. ] 

CIO]  A 

Cll]  A   Y. KATZIR. I. A.F. .NOVEMBER  1987 

C12]  A 


PC«-10000 
FS«-50000 

TC«-FS  DIGFREQ  FC 


TC 
X.  256637061 

TCi-PI 
0.4 
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7LPCOEFFC03V 

CO]  hLP«-n  LPCOEFF  TC:I;0IO 

C13  A 

C2J  A   CALCULATING  THE  'n'  LOWPASS  FIR  FILTERS  COEFFICIENTS  USING  FOURIER 

C3  3  A   METHOD 

C4:  a 

C5i  aio«-i 

C6J  n«-n-l 

C71  n«-Un+2) 

C8J  hLPHipO 

C9]  hLP^(l+(PIxCinn)xSIN(Cl.n)xTC) 

C10J  hLP«-(TC*PI).hLP 

£11]  hLP«-(nTC*hLP)).hLP 

C12  3  QIO«-0 

C131  -*0 

C143  A 

C15  3  A   THE  COEFFICIENTS  ARE  FOR  CAUSAL  NONRECURSIVE  FILTER. 

CIS]  A   THE  COEFFICIENTS  ARE  NORMALIZED. 

C17]  A   n«-ODD  NUMBER  OF  REQUIRED  COEFFICIENTS  C RECTANGULAR  WINDOW) 

C13:  A   TC«-DIGITAL  CUTOFF  FREQUENCY  CRAD.  2 

C191  A 

C20J  A   Y. KAT2IR. I. A. F. .NOVEMBER  1987 

C211  A 


TO-0.4XPI 
n«-21 


hLP*n  LPCOEFF  TC 

hLP 
"1.559907461E~17  "0.0336367435  "0.02338723209  0. 0267282652S  0.05045511524 
"1.559907461E~17  "0.07568267286  "0.06236595225  0.09354892838 
0.3027306915  0.4  0.3027306915  0.09354892838  "0.06236595225 
"0.07568267286  "1. 559907461E"17  0.05045511524  0.02672826525 
"0.02338723209  "0.0336367435  "1.559907461E"17 


76 


▼HPCOEFFCCn? 

COJ  hHP«-n  HPCOEFF  TC:QIO:N:hLP 

CIJ  a 

C2J  A   CALCULATING  THE  *n'  HIGHPASS  FIR  FILTER  COEFFICIENTS  USING  FOURIER 

C3]  A   METHOD 

C-U  A 

cs]  aio^o 

C61  N«-n 

C7]  TO-PI-CTC3 

C8]  hLP+n   LPCOEFF  TC 

C93  hHP«-CC-l)*l.N)xhLP 

CIO]  -»0 

cn:  a 

C12]  A   THE  COEFFICIENTS  ARE  FOR  CAUSAL  NONRECURSIVE  FILTER. 

C13]  A   THE  MAGNITUDE  IS  NORMALIZED. 

C14]  A   n«-ODD  NUMBER  OF  REQUIRED  COEFFICIENTS  (RECTANGULAR  WINDOW) 

C15]  A   TC*-DIGITAL  CUTOFF  FREQUENCY  CRAD. ] 

CIS]  A 

C17]  A   Y. KATZIR. I. A.F. .NOVEMBER  1987 

CIS!  A 


TCWJ.4XPI 
n«-21 

hHP«-n  HPCOEFF  TC 

hHP 
■2.339861191E~17  0.0336367435  0.02338723209  "0.02672826525  "0.05045511524 
-2.339861191E"17  0.07568267286  0.06236595225  "0.09354892838 
"0.3027306915  0.6  "0. 3027306915  "0.09354892838  0.06236595225 
0.07568267286  "2. 339861191E"17  "0.05045511524  "0.02672826525 
0.0233S723209  0.0336367435  "2.339861191E"17 
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vbpcoeffcchv 

CO]  hBP«-n  BPCOEFF  TUL ; OIO ;TC; TO ;hLP 

CI]  A 

C2]  A   CALCULATING  THE  'n'  BANDPASS  FIR  FILTER  COEFFICIENTS  USING  FOURIER 

C3]  A   METHOD 

C4J  a 

C5]    0IO«-0 

C63    TC«-<TULC0]-TULCl]]+2 

C7]    TO«-CTULC0]*TULCl])+2 

C81    hLP«-n  LPCOEFF  TC 

C9]    n«-n-l 

CIO]   n»l(n+2) 

Cll]   HLP*-nihLP 

C12]      HBP«-C2xCOSCUn*l]xTO))xhLP 

C13]  hBP«-CnfC«hBP)).&BP 

C14]  -*0 

CIS]  A 

C16]  A   THE  COEFFICIENTS  ARE  FOR  CAUSAL  NONRECURSIVE  FILTER. 

C17]  A   THE  COEFFICIENTS  ARE  NORMALIZED. 

C18]  A   n«-ODD  NUMBER  OF  REQUIRED  COEFFICIENTS  (RECTANGULAR  WINDOW] 

C19]  A   TUL«-DIGITAL  UPPER  AND  LOWER  CUTOFF  FREQUENCIES  CTWO  ELEMENTS)  CRADJ 

C20]  A 

C21]  A   Y. KATZIR. I. A. F. .NOVEMBER  1987 

C22]  A 


TUL«-2p0 
TULCl]«-0.6xpI 
TULC2]«-0.4xPI 
W21 

hBP<-n   BPCOEFF  TUL 

hBP 
.799537305E-18    1. 204383817E~17    0.04677446419    "3.152585957E"17    "0.1009102 
3.39562730SE"17    0.1513653457    "3. 156340459E-17    "0.1870978568 
1.207227699E"17    0.2    1.207227699E~17    "0.1370978568    "3. 156840459E"17 
0.15136534S7    3. 395627305E~17    "0.1009102305    "3. 152585957E~17 
0.04677446419    1. 204383817E"17    ~7.799537305E~18 
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7BSCOEFFC037 

CO]  hBS«-n  BSCOEFF  TUL;aiO:hBP 

CI]  a 

C2  3  A   CALCULATING  THE  'n'  BANDSTOP  FIR  FILTER  COEFFICIENTS  USING  FOURIER 

C3  3  A   METHOD 

C43  A 

C5  3  aio«-o 

CSI  hBP«-n  BPCOEFF  TUL 

C73  hBS«-hBP 

C8I  n«-n-l 

C93  hBSCn+23«-l-hBPCn+2] 

C103  -*0 

1113  A   THE  COEFFICIENTS  ARE  FOR  CAUSAL  NONRECURSIVE  FILTER. 

C12  3  A   THE  COEFFICIENTS  ARE  NORMALIZED. 

C13  3  A   n^ODD  NUMBER  OF  REQUIRED  COEFFICIENTS  (RECTANGULAR  WINDOW) 

C14  3  A   TUL«-DIGITAL  UPPER  AND  LOWER  CUTOFF  FREQUENCIES  CTWO  ELEMENTS)  CRADJ 

C153  A 

C163  A      Y. KATZIR. I. A.F. .NOVEMBER   1987 

C173  A 


TULC13«-0.4xPI 
TULC23*-0.2xPI 
n«-21 

hBS«-n   BSCOEFF   TUL 

UBS 

7.799537305E"18    0.01284809274    "0.01445410434    "0.06997550689    "0.08163809137 
2.341241641E"17    0.1224571371    0.1632761827    0.05781641735    "0.115632834 
0.3    "0.1156328347    0.05781641735    0.1632761327    0.1224571371 
2.341241641E"17    "0.08163809137    "0.06997550639    "0.01445410434 
0.01284809274    7.799537305E~18 
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▼HAMMINGCtm 

CO]  W«-HAMMING   M:M;3I0 

CI]  A 

C2]  A      GENERATING  THE    'N'    WIDE   HAMMING  WINDOW 

C31  A 

C4i  aio«-o 

C5]  W«-HpO 

C6]  -»((N+2]^(rCN+2))J/ODD 

C7]  A  EVEN  HUMBER  OF  SAMPLES  IN  WINDOW 

C8]  M«-CN+2) 

C9]  W«-0.54t-0.46xCOSCCCPIx2)x\.M)+W5 

CIO]  W«-(eW).W 

Cll]  *0 

CI2]  0DD:M*-(N*l)+2 

C13]  W«-0.54*0.46xCOSCCCPIx2)xUi]+HJ 

C14]  W«-<eW).UW 

C15]  -*0 

C16]  A 

C17]  A   THE  WINDOW  STARTS  AT  t=0  AND  HAS  'N'  SAMPLES. 

C18]  A   N+NUMBER  OF  SAMPLES  C VECTOR] 

C19]  A 

C20]  A   Y. KATZIR. I. A. F. .NOVEMBER  1987 

C21]  A 


N«-31 

W4-HAMMING  N 

W 
0.08236011124  0.1010959471  O-1378005566  0.1909712636  0.258431248  0.3374186 

0.4246998352  0.5167013323  0.6096567777  0.6997604163  0.7833234448 

0.8569247828  0.917551183  0.9627205933  0.990583773  1  0.990583773 

0.9627205933  0.917551183  0.8569247828  0.7833234448  0.6997604163 

0.6096567777  0.5167013823  0.4246998352  0.3374186903  0.258431248 

0.1909712636  0.1378005566  0.1010959421  0.08236011124 
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7HANNINGCn]V 

CO]  W«-HANNING   N:M:0IO 

CI]  A 

C2]  A   GENERATING  THE    'N'    WIDE   HANNING  WINDOW 

:n  a 

C4J      aio«-o 

C5  3    W«-Np0 

C6]    -KCN+2)*CrCN+2n)/ODD 

C71  A  EVEN  NUMBER  OF  SAMPLES  IN  WINDOW 

C8]    M«-CN+2) 

C91    W«-0.5*0.5xCOSCCCPIx2)xi.M)+N) 

CIO]   W«-(eW).W 

C113   -»0 

C12]  ODD:M«-CN+l]+2 

C13]   W*O.5+0.5xCOSCCCPIx2)xi.M)+N) 

C14]   W«-CeW).liW 

CIS]   -»0 

C16]  A 

C17]  A   THE  WINDOW  STARTS  AT  t*0  AND  HAS  'N'  SAMPLES. 

C18]  A   N«-NUMBER  OF  SAMPLES  C VECTOR] 

C19]  A 

C20]  A   Y. KATZIR. I. A- F. .NOVEMBER  1987 

C21]  A 


N«-31 


W4-HANNING  N 


W 


2.565338304E~3  0.0229303718  0.06282669193  0.1206209387  0.1939470087 


0.2798029242  0.3746737339 
0.76448200S2  0.3444834595 
0.9897649706  0.9594789058 
0.6736526264  0.5757138888 


0.4746754156  0.5757138888  0.6736526264 
0.9103817206  0.9594789058  0.9897649706 
0.9103817206  0.8444834595  0.7644820052 
0.4746754156  0.3746737339  0.2798029242 


0.1939470087  0.1206209387  0.06282669193  0.0229303718  2. 565338304E"3 
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VBARTLETTCaJV 

£0]  W«-BARTLETT  N:M:£JIO 

Cll  A 

£2]  A   GENERATING  THE  'N'  WIDE  BARTLETT  (TRIANGULAR)  WINDOW 

C3]  A 

C4]  aio*o 

£5]  W«-NpO 

£6]  -KCN+2D*CrCN+2)))/ODD 

£7]  A  EVEN  NUMBER  OF  SAMPLES  IN  WINDOW 

£8]  M«-<N+2) 

£93  W«-(e(l-(Ui)+Mn.l-CU«+M 

£10  3  -»0 

£113  0DD:M«-CN*l)+2 

£121  W«-1-(\.M)+M 

£13  3  W«-(eW),UW 

£14]  -*0 

£15]  A 

£16]  A   THE  WINDOW  STARTS  AT  t»0  AND  HAS  'N'  SAMPLES. 

£171  A   N«-NUMBER  OP  SAMPLES  C VECTOR] 

£13]  A 

£191  A   Y. KATZIR. I. A. F. .NOVEMBER  1987 

£20]  A 


N«-31 

W4-BARTLETT  N 

W 
0.0625  0.125  0.1875  0.25  0.3125  0.375  0.4375  0.5  0.5625  0.625  0.6875  0.75 

0.8125  0.375  0.9375  1  0.9375  0.375  0.3125  0.75  0.6875  0.S25  0.5625  CT$ 
0.4375  0.375  0.3125  0.25  0.1875  0.125  0.0625 
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VFREQRESCD]* 

CO]  HW«-M  FREQRES  h:OIO:N:H 

ClJ  a 

C2]  A   GENERATING  FREQUENCY  RESPONSE  MAGNITUDE  AND  PHASE  FROM  TIME  DOMAIN 

C3]  A   COEFFICIENTS 

C41  A 

ts]   aio«-o 

C6J    N«-CM-ph)pO 

C7J    h«-h.N 

C8]    H«-FFT  h 

191        HW«-(M+2)TCXMAGN  H) 

CIO]   PW«-<M+2)TCXPHAS  H) 

til]   -»0 

C12]  A   THE  FUNCTIONS  FFT.XMAGN.XPHASE  ARE  USED. 

C13]  A   IN  ORDER  TO  RECEIVE  A  BETTER  RESOLUTION  THE  COEFFICIENTS  ARE  2ERO- 

C14]  A   PADDET  TO  'M'  SAMPLES. 

CIS]  A   M«-NUMBER  OF  SAMPLES  INCLUDING  ZERO  PADDING  CUSE  ONLY  RADIX  TWO  NO. 

C16]  A   h*-THE  FILTER  COEFFICIENTS  C VECTOR] 

C17]  A   HW: FREQUENCY  RESPONSE  MAGNITUDE  CM+2  SAMPLES) 

CIS]  A   PW: FREQUENCY  RESPONSE  PHASE  (M+2  SAMPLES)  CRAD. ] 

C19]  A 

C20]  A   Y.KATZIR. I. A. F.NOVEMBER  1987 

C21]  A 


HLP*21   LPCOEFF    C0.4xpi) 

HLP*-128    FREQRES   hLP 

HLP 
0.9567807992  0.9619585014  0.9763182673  0.9965906571  1.018111173  1.03582952 
1.045421125  1.044271236  1.032124069  1.01125183  0.986091438  0.9624009 
0.9460845406  0.9419036347  0.9523233652  0.9767236974  1.011145443 
1.048645824  1.080227975  1.096202394  1.087755987  1.048461283  0.975462 
0.3701304604  0.7380554233  0.5883839646  0.4325882327  0.232876649 
0.1505027492  0.04424706035  0.03068671843  0.0731944037  0.08614789716 
0.07567134118  0.04996228114  0.01787457278  0.01248921845  0.0349417940 
0.04595163274  0.04486445606  0.03360579506  0.01595981635  3 . 413311717E"3 
"0. 02001949514  0.03043376941  0.03290460379  0.0275668884  0.01625216492 
1.963679943E~3  0.01185437487  0.02209076454  0.02662206686  0.024719091 
0.01713279563  5. 357126075E"3  6. 366569896E"3  0.01668743018  0.02278343 
0.02336291733  0.01842504363  9.228990525E_3  2. 014864013E"3  0.01265918 
0.02022375085 
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VIDEALFCQ]? 

[01  H*-TYPE  IDEALF  TC;DIO 

ci]  a 

C2J  R   GENERATING  IDEAL  FREQUENCY  RESPONSE 

[3J  A 

C4J  0104-0 

[5]  H*-TC[0]pO 

[6]  -KCx/CTYPE='LP'nv(x/CTYPE»'HP')n/LP 

[7]  -»<(x/(TYPE»'BP' ))vCx/CTYPE='BS' )))/BP 

[8]    'IMPROPER  ARGUMENT. USE  ONLY:  LP.  HP  .  BP  .OR   BS  . 

C9I  -»0 

CIO 3  LP:HCULCCTCC0]+2)xTCtl])*l]]«-l 

Cll]  H«-(TC[0]+2)TH 

[12]  -»(x/ (TYPE* 'LP')) /CONTINUE 

[13]  H«-~H 

[14]  -*C0NTCNUE 

115 J  BP:H[CUUCTCC0]+2)x(TC[l]-TC[2]  ))))♦(!.  CCCTCC0]+2]xTCC2])*m]«-l 

[16]  H«-(TC[0]+2)TH 

[17]  -»<x/CTYPE='BP' )) /CONTINUE 

[18]  H«~H 

[19]  CONTINUE: H«-H.*H 

[20]  -»0 

[21]  A 

C22J  A   'TYPE'«-THE  TYPE  OP  FILTER  TO  BE  GENERATED  CLP.HP.BP.BS J 

[23]  A   TCCO]«-NUMBER  OF  SAMPLES  CUSE  ONLY  RADIX  TWO  NO.) 

C24]  A   TC[1]«-THE  CUTOFF  DIGITAL  FREQUENCY  FOR  LP  OR  HP  FILTER  OR 

[25]  A   THE  UPPER  DIGITAL  FREQUENCY  FOR  BP  OR  BS  FILTERS 

[26]  A   TC[2]«-THE  LOWER  DIGITAL  FREQUENCY  FOR  BP  OR  BS  FILTERS 

[27]  A   THE  FREQUENCIES  SHOULD  BE  A  FRACTION  OF  PI  [RAD.]. 

[28]  A 

[29]  A   Y.KATZIR  .I.A.F. .NOVEMBER  1987 

C30]  A 
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0   0 


TCtl3«-128 
TCC2J«-0.6 
H«-'HP'    IDEALF  TC 

a 

0000000000 

liiiiiiii 
111111111 

000000000 

0 

1 
1 

0 

0 

1 
1 

0 

0 

1 
1 

0 

0 

1 

0 
0 

0 

1 

0 
0 

0 

1 

0 

0 

1 

0 

0 

1 

0 

0    0    0    0    0    0    0 

1111111 

0    0    0    0    0    0    0 

0    0    0    0    0 

11111 

0    0    0    0    0 

0    0    0    0 

1111 

0    0    0    0 

0 

1 

0 

TO3p0 

TCC1J«-128 
TC£23«-0.3 
TCC3J«-0.12 


H«-'BP'    IDEALF  TC 


H 

000000001111111111100 
000000000000000000 
000000000000000000 
11111100000000 


0000000000000000 
0000000000000000 
0000000000000011 
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VFCOEETtrm 

CO]  h«-n  FCOEFF  HW;OIO:H:m;SIZE 

:i:  a 

C2]  A   CALCULATING  THE  'n'  FILTER  COEFFICIENTS   USING  IFFT 

C3]  A 

C4]    010*0 

C5J   n*n-l 

C6]    SIZE«-Un+2J 

C71   nn-SIZE+1 

C81    H>IFFT  HW 

C93    h«-(l.m)TH 

CIO]  h«-mph 

Cll]   h*-CSIZETC<J>h)).h 

C12]   -»0 

C13]  A 

C14]  A   THE  COEFFICIENTS  ARE  FOR  CAUSAL  NONRECURSIVE  FILTER. 

C15  3  A   n«-ODD  NUMBER  OF  REQUIRED  COEFFICIENTS  C RECTANGULAR  WINDOW) 

C16]  A   HW«-  THE  IDEAL  FILTER  SAMPLES  (RADIX  TWO) 

C171  A 

C18]  A   Y.KATZIR.I.A.F. .NOVEMBER  1987 

C19]  A 


0  0 


H 
0000001111111111 
000000000000000 
000000000000000 
11111100000000 

10  0  0  0  0  0 
0  0  0  0  0  0  0 
0  0  0  0  0  0  0 

00000000 
00000000 
00000000 

0  0  0  0 
0  0  0  0 

0  0  11 

n*-2l 

h«-n  FCOEFF  H 

h 

0.02511120991  0.04315473177  0.03628637547  "6.290075645E~3  "0.07065358471 

"0.1217862264  "0.1233903877  "0.06216894475  0.03973030337  0.133830098 
0.171875  0.1338300989  0.03973030337  "0.06216894475  "0.1233903877 
"0.1217862264  "0.07065358471  "6.290075645E"3  0.03628637547  0.0431547 
0.02511120991 
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)LOAD  0  EC3410 
0  EC3410  SAVED  12/11/1987  12:48:01 


NAVAL  POSTGRADUATE  SCHOOL 


MONTEREY . CALIFORNIA 


ECE  -3410  WORKSPACE 


THESIS  BY 
Y.KATZIR.I.A.F. 
ADVISOR: PROF.  C.W.  THERRIEN. 


VERSION  1.0 


SEPTEMBER  1987 


)FNS 

CC 

COS 

COV 

FCC 

FFT 

FLCV 

IFFT 

LPFP 

MEAN 

PI 

PSE 

PSEB 

SACF 

SIN 

XPHAS 

XSQRT 

XTIME 

INFORMATION      LCV 
XCONJ    XEXPO    XMAGW 
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VSACFCD3V 
CO]      R«-LAG   SACF    DAT ;DIO: MEAN ;N 

[2  3       «    SAMPLE   AUTOCORRELATION   FUNCTION    CUNBIASED) 

C33  A 

[4]  DIO*-K«-0 

[53  N+-pDAT 

[63  R<-LAGp0 

[7  3  DAT«-DAT,R 

C83  L00P:R[K3*-O/DATxC-KD<J>DAT)*N 

[93  -KLAG>K«-K+l)/LOOP 

CIO  3  -»0 

C123    A    Y.KATZIR.I.A.F..JULY    1987 


[133  A 

C143  A 

C15  3  A 

C163  A 


C14  3    A    DAT«-DATA  VECTOR    CONE    DIMENSIONAL) 
C15  3    A    LAG«-NO.    OF    LAGS    (MAX.    3  00) 


RXX*-100    SACF   GAUSD 


Rxx 

1.1816    "0.011604    "0.066841    "0.14886    "0.11086    "0.0S0135    "0.28244 

0.12983    0.011973    *2.0869E~3    "0.052813    "0.066688    0.1414    0.16079 
0.021921    "0.050052    "0.03314    "0.12319    "0.11111    "0.087361 
7.9253E"3    0.047914    0.082207    0.022012    "0.072971    0.047907 
"0.028347    0.03337    "0.091018    0.11214    4.3453E"3    "0.020604    0.08519 
"0.13963    0.038409    "0.059126    "0.023359    "0.04759    0.03288    0.23321 
"0.021454    0.051991    "0.049765    "0.025113    0.059572    "0.11704    0.11551 
"0.11528    "0.058985    "0.038291    "0.077532    0.11274    0.10221    0.18426 
"0.09541    "0.14171    8.3982E"3    "0.10611    "0.13464    0.037944    0.035897 
0.13416    "0.081343    "0.028463    0.13128    0.026911    0.011188    "0.048038 
0.048817    "0.023835    "0.14964    "0.071093    0.084411    0.04223    0.055516 
0.08221    "0.011808    "0.03549    "4.7687E"3    "0.011637    0.010587 
*7.3499E"3    0.033257    "0.065615    "0.036728    0.038447    0.014501 
0.029666    0.042598    "0.096824    5.4581E"3    0.052421    0.047676 
"0.044033    "0.049402    0.016731    "1.8189E"3    6.5336E"3    "0.028012 
0.015022 
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SO  MEAN  UHIT2  SAUSSIriN  N0IS2 


3. 

2. 

I. 

'I 
A 
II     0. 


I A  in 

Pi  it   i  I A 

ffrwj  ii  r 


L  !i 


'M-'ihi-l.Triiiii'ynlrvr  l-riu/i  ■•? 


Mii!-a 


!j  3 


20  40  60 

TIKE  3WL2S 
ZZ2G  MEAN  1*172  UNIFORM  NOISS 


30 


ICC 


hi 


I 

I 
I  y 

.  i    l 

f-fc-  - 1 — If 


f         .     I 
I 


I VI 


I 


i  i  •  hi 


V 


I 


I: I  I   RIH 


Mm*-- 

*f  [M-   •r.»-   -•'    i-(  v»— i 


IE 
1 

!!  H 


20 


40  SO 

TIME  SAMPLZ3 


30 


fl 


100 


COBKIATION  OF  >*ITZ  H0I32. 


0.6 


'0.3  t- 


-i 

i 

-i 


o  %  k  i\^^^ /I^Y^-V^ 


-0.3 


111  ..I 


20 


40  60 

lac; 


30 


100 
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VMEANCOJV 

CO]  M«-MEAN  F 

CI  J  fl  CALCULATE  THE  MEAN 

C23  M<-(+/  +  /F)tx/pp 

C3]  n 

£4  3  b  SOURCE :MDSP  W. SPACE 

C53  n  F*-DATA  VECTOR 

C63  a  M:THE  MEAN  (NUMBER) 

C7]  r 


X«-tlOO 

X 

12  3  4  5  6 


248529  30%!  g  S  "  J2  JJ  375  JJ  39  2  X  J!  "  "  "  24  »  "  27 

50  51  52  53  54  55  56  57  58  59  60  11    «  2  2  2  44  45  46  47  48  49 

72  73  74  75  76  77  78  79  80  81  M  2  2  2   4  "  M  6?  68  69  70  71 

94  95  96  97  98  99  100  "  84  M  86  8?  88  89  90  91  92  93 


50.5 


XX<-MEAN  X 
XX 

x<-x-xx 

X 


"49.5    "48.5    -47.5      46.5    "45.5    "44.5    "43    5    -42    e    -,, 
-JS.5    "37.5    "36.5    "35.5    "34    5    41   5    "32*1    -2 

fti  ft!  ft!  33  %■  tl  fl  f-|  a..  3.:. 


"7.5    -6.5    -5.5    -4.5    -3    5    "2   5  -t  V  -n«*   ^^    ~10'5    ~9.5    "8.5 
6-5    7.5    8.5    9.5    lofs    li   5    12   5    13    s    il  V^    c*5    2*5    3*5    4"5    5"5 

19.5    20.5    21.5    22.5    23    5    24    5  25    5 %J    ^    ,^;5    16'5    17*5    18-5 

32.5    33.5    34.5    35.5    36    5    37*1  U J    J?'!    f Z * f    28'5    29-5    30.5    31.5 

45.5    46.5    47.5    48.5    49.5  39*5   4°'5    41-S    42'5    43.5    44.5 

XX«-MEAN  X 
XX 
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VCOVCD]? 

C03  k«-cov  R;aio;N:i 

ci]  a 

[23  A   GENERATE  TOEPLIT2  COVARIANCE  MATRIX  FROM  COVARIANCE  FUNCTION 

[33  a 

[43    CIO*-3>0 

[53    N«-pR 

[63    K«-CN.N3pO 

[73    R*-(4>liR3.R 

[83  LOOP:KCI:3*C-N)TC-I3AR 

[93    -»CN>I«-I  +  l)/LOOP 

[103  A 

[113  A   SOURCE:  MDSP  W. SPACE 

[123  A   R<-THE  CORRELATION  FUNCTION  (VECTOR) 

[13  3  A 


RCLl*-3  SACF  CL1 


RCL1 
7.3328  6-8337  6.4035 


CMCLK-COV  RCL1 


GMCL1 
7.3328  6.8337  6.4035 
6.8337  7.3328  6.8337 
6.4035  6.8337  7.3328 
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?LPFPCD]V 

CO]   a<-P  LPFP  DATA;DIO:R:K:MDATA;M 
CI]   A 


C2]  fl  GENERATE  THE  Pth  ORDER  LINEAR  PREDICTION  FILTER  PARAMETFRe; 

C3]  fl  AND  THE  PREDICTION  ERROR  VARIANCE  PARAMETERS 

C4]  A 

C5]    DIO+-0 

C6]    a«-CP+13pO 

C7]    M«-MEAN  DATA 

C8]  MDATA«-DATA-M 

C9]  R<-CP+13SACF  MDATA 

CIO]  K«-COV  R 

Cll]  a«-i.a«~  CliR3BCP,P)TK 

C12]  PERR«-+/CaxR) 

C13]  -»0 

C14]  fl 

III]    fl  F^cSS"0N  S°LVES  N0RMAL  EQUATI0NS  USING  THE  SACF.COV.AND  MEAN 
C17]  fl 

C18]  fl  Y. KAT2IR. I. A.F. .AUGUST  87 
C19]  fl 

C20]  fl  a:  FILTER  COEFFICIENTS 

\ll\  "  !E^  THE  ESTIMATE  PREDICTION  ERROR  VARIANCE 

C22]  fl  P<-THE  FILTER  ORDER 

C23]  fl  DATA«-THE  DATA  VECTOR 

C24]  fl 


a<-3  LPFP  CL1 


a 
1  "0.897346235  "0.017308669  "0.02113346634 


PERR 
0.9625826774 
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C03  Y«-X  CC  H:DIO:N;I:NN;SNN 

CI]  ft 

C2]  ft    1-D   CIRCULAR   CONVOLUTION 

[33  ft 

C4]         DIO«-I<-0 

C5]  SNN*-pCNN«-pH) 

C6]  N<-pX 

C71  -KN*NND /ERROR 

C8]  -KSNN*  13 /ERROR 

C9]  Y<-NpO 

CIO]  LOOP:YCI]«-«-/XxC-l+IDO<t>H 

Cll]  ■+CN>I<-3>l)/LOOP 

C12]  -»0 

C13]  ERROR: 'ARGUMENTS  ARE  NOT  THE  SAME  SIZE  OR  1-D. 

C14]  -»0 

C15]  ft 

C16]  A  Y.XATZIR. I. A. F.AUGUST  1987 

C17]  ft 

C18]  ft  ARGUMENTS  NEED  TO  BE  THE  SAME  SIZE  AND  1-D. 

C19]  ft  X  AND  H<-ARGUMENTS  TO  BE  CONVOLVED  CVECTOR] 

C20]  ft 


3  2  1 
3  2  1 


X 


H 


Y<-X  CC  H 

Y 

13  13  10 

LONG 
1.23  '4  56789  10 

X'CC  LONG 
ARGUMENTS  ARE  NOT  THE  SAME  SIZE  OR  1-D. 
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vFcctmv 

CO]   Y<-X  FCC  H;aiO;N;NN;SNN;A 

CI]   A 

C2]   A  FAST  1-D  CIRCULAR  CONVOLUTION 

C3]   A 

[4]    0IO<-0 

C5  3    SNN<-pCNN<-pHD 

[6]    N«-pX 

[7]         -KNN>80)/STATM 

[8]         -KN*NN]/ERROR 

[93    -KSNN*1)  /ERROR 

CIO]   A«-CN.N3p«H 

C113      A<-C-vNDOA 

C12]      Y<-l<J>CA+.xX) 

[13  3      -»0 

[143  ERROR:  'ARGUMENT  ARE  NOT  THE  SAME  SIZE  OR  1-D. 

[153  -+0 

C16  3  STATM: 'PLEASE  USE  CC  OR  LCV  FUNCTIONS.' 

C173   -»0 

C133  A 

C19]  A  Y. KATZ IR.T. A. F. .AUGUST  1987 

C20]  A 

C213  A  ARGUMENTS  NEED  TO  BE  THE  SAME  SIZE  AND  1-D. 

C22  3  A  LENGTH  OF  ARGUMENT  VECTOR  NEED  TO  BE  SHORT. 

[23]  A  X  AND  H«-ARGUMENT  TO  BE  CONVOLVED 

[24]  A 


X 

3  2  1 

H 
3  2  1 

Z«-X  FCC  H 

Z 
13  13  10 
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7Lcvcn]v 

£03  Y«-X1  LCV  Hl:DIO;Nl:N2;N;XZ:HZ;X;H 

CI]  A 

[2]  A  1-D  LINEAR  CONVOLUTION 

[3]  A 

[4]    QIO<-0 

[53  Nl<-pXl 

[6]  N2<-pHl 

[7]  N<-N1+N2-1 

[8]  XZ«-NTX1 

[9]  HZ«-NTH1 

CIO]  Y«-XZ  CC  HZ 

[11]  -*0 

1121  A 

C13]  A  LINEAR  CONVOLUTION  USING  CC  FUNCTION 

C14]  A 

C15]  A  Y.KATZIR.I.A.F. .AUGUST  1987 

C16]  A 

C17]  A  X  AND  H<-THE  ARGUMENTS  TO  BE  CONVOLVED  C VECTORS] 

C18]  A 


3  2  1 
3  2  1 


X 


H 


Y«-X  LCV  H 

Y 

9  12  10  4  1 

LONG 

1  2  3  4  5  6  7  8  9  10  11  12  13  14  15  16  17  18  19  20  21  22  23  24  25  26  27  28  19 

30-31  32  33  34  35  36  37  38  39  40  41  42  43  44  45  46  47  48  49  50  51  52  53 

54  55  56  57  58  59  60  61  62  63  64  65  66  67  68  69  70  71  72  73  74  75  76  77 

78  79  80  81  82  83  84  85  86  87  88  89  90 

Y«-X  LCV  LONG 

Y 

3  8  14  20  26  32  38  44  50  56  62  68  74  80  86  92  98  104  110  116  122  128  134  140 

146  152  158  164  170  176  182  188  194  200  206  212  218  224  230  236  242  248 

254  260  266  272  278  284  290  296  302  308  314  320  326  332  338  344  350  3S"4 

362  368  374  380  386  392  398  404  410  416  422  428  434  440  446  452  458  464 

470  476  482  488  494  500  506  512  513  524  530  536  269  90 
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VFLCVCD3V 

CO]  Y<-X1  FLCV  H1;DI0;N1;N2:N;XZ:HZ:X:H 

C1J  a 

121  A  FAST  1-D  LINEAR  CONVOLUTION 

[3]  A 
C4  J    DIO«-0 
[5]    Nl«-pXl 
[6]    N2<-pHl 
C7]    N«-N1+N2-1 
[8]    XZ*-PTX1 
[9]    HZ«-NTH1 
CIO]   Y<-XZ  FCC  HZ 
Cll]   -»0 

C121  A 

[13]  A  FAST  LINEAR  CONVOLUTION  USING  FCC  FUNCTION 

[14]  A 

[15]  A  Y.KATZIR.I.A.F- .AUGUST  1987 

[16]  A 

[17]  A  X  AND  H«-THE  ARGUMENTS  TO  BE  CONVOLVED  C VECTORS) 

[IS]  A 


X 

3  2  1 

H 
3  2  1 

Y«-X  FLCV  H 

Y 

9  12  10  4  1 


Y«-X  FLCV  LONG 
PLEASE  USE  CC  OR  LCV  FUNCTIONS. 
VALUE  ERROR 
FLCV[10]  Y<-XZ  FCC  HZ 
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VFFTCCHV 
CO]   z«-FFT  x;k:n;r;QIO 
CI]    QICH-0 

C2]    -Kr*-C2€ppX)A2€lTpX)T8 
C3]    -Kl€pk)T7-r«-25>k«-L0.5xpx 
C4]    'NOT  A  COMPLEX  VECTOR' 
C53    -*0 

C6]    -»8,px«-x.[~0.5]0 
C7]    x«-$Ck.2)pX*2 

C8]    -K20ak«-  "1  2  4  8  16  32  64  128  256  512  1024  2048  4096  8192  16384  3276 
C9]    -KHpZ«-xU  0  0 
CIO]   'NOT  A  POWER-OF-2  ORDER' 
Cll]   -»-0EX  'Z' 
C12]   z«-Cl  2  ,kp2)px 
C13]   -fki  0  28.23 
C14]   -»(2€0NC  'wt'U16 
C15]   -Kk€ppx«-wt)T20 
C16]   n«-0.-<&Hx«-2o(O2*n)x\.n*4 
C17]   wt«-x«-Ckp2)p(x,n)  ,n.-x 
C18]   -»20 

C 19  ]   x«-« (kp  '  xC  ;;;;;:;::;;;;;:;;:;').'  0  ] ' 
C20]   n*— /[2)z 
C21]   k<-pz«--t-/[2]z 

[22]   -K4<k*-ppz«-Z.[0]C-/[l]nxkpx).[0.5]  +  /[l]nxkpex)T19 
[23]   n«— /[2]z 
[24]   Z«-+/[2]Z 
[25]  nC:0:lJ«~nC;0:l] 
[26]   n[;  1  0  ;l]«-n[;:l] 
[27]   z«-z,[0)n 
[28]   z«-*C+/z),[n*-0]-/z 
[29]   -n.r 
[30]  n*-~!Q<t>z 
[31]   n[l;]«--n[l;] 
[32]   x<-n-z 
[33]   n«-~ltpz«-z+n 
[34]   k<-0,<t>lir«-2o(o-rn)x\.n*2 
[35]   r«-.(pz)pk,r,-r 
[36]   x«-(+j*xxr)  ,["0.5]-/rxex 
[37]   z«-Cz.+x)  ,z-x 
[38]   -»0 
[39]  A 

[40]  A  APL-STYLE  FAST  FOURIER  TRANSFORM 

[41]  fl   X«-CAN  BE  EITHER  REAL  CA  VECTOR)  OR  COMPLEX  C2-ROW  MATRIX) 
[42]  fl   OPTIMIZED  FOR  SPEED 

[43]  B   Z:THE  OUTPUT  IS  COMPLEX  (REAL  AND  IMAGE  PARTS). 
[44]  A   SOURCE: PAUL  PENFIELD, JR. 6/13/79 
[45]  ft 
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VIFFTCQ]? 

CO]  z<-IFFT  x;k;n;r;DIO 

ci]   nio*-o 

C2]         -Kr«-C2€ppx]A2€lTpx]T8 

C3]         -Kl€pk]f7-r<-25>k<-l0.5xpx 

C4]  'NOT   A   COMPLEX   VECTOR' 

C53         -»0 

C6]         -»8,px«-X.C  "0.510 

C7]         x«-«Ck.2)px*4 

C8]    -»C20£k«-  "1  2  4  8  16  32  64  128  256  512  1024  2048  4096  8192  16384  3276 

C9]    -Klipz«-xU  0  0 

CIO]   'NOT  A  POWER-OF-2  ORDER' 

Cll]   •♦-DEX  '2' 

C12]   z«-Cl  2  ,kp2]px*n 

C13]   -»ki  0  28  23 

C14]   -»(2€0NC  'wt'  U16 

C15]   -Kkeppx«-vt)T20 

C16]   n«-0.-<t>lAx<-2O(O2-Ml)xT.n-5-4 

C17]  vt«-x«-ckp2]p(x.n]  ,n,-x 

C18]   -+20 

C19]  x*-«ckp'xC: :;::::;::::::::::;').  '0] ' 

C20]   n«—  /C2]z 

C21]   k*-pz<-+/C2]z 

C22]   -»(4<k*-ppz«-z.C0]C  +  /Cl]nxkpx).C0.5]-/Cl]nxkpe-x]T19 

C23]   n«--  /C2]Z 

C24]   Z«-+/C2]z 

C25]  nC:l:l]<— nC:i:i] 

C26]   nC;  1  0  ;13«-nC:;l] 

C27]   z«-z,C0]n 

C28]   z<-«(  +  /z]  ,Cn«-0]-/z 

C29]   -»ir 

C30]      n<-_l<t><t>z 

C31]   nCl;]«--nCl;] 

C32]   x«-z-n 

C33]   n<--lTpz«-z+n 

C34]   k*-0,<01ir«-2oCo-rn]xT.n-i-2 

C35]   r«-(pz)pk,r.-r 

C36]   xt-C-j'xxr]  ,C"0.5]+/xxer 

C37]   z«-(z+xJ,z-x 

C38]  -»0 

C39]  fl 

C40]  A   APL-STYLE  INVERSE  FOURIER  TRANSFORM 

C41]  ft   X<-CAN  BE  EITHER  REAL  CA  VECTOR]  OR  COMPLEX  C2-ROW  MATRIX] 

C42]  A   OPTIMIZED  FOR  SPEED 

C43]  A   SOURCE: PAUL  PENFIELD.JR.  6/13/79 

C44]  ft 
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X 

1111111111111111 


Y<-FFT   X 

Y 

16 

0       0       0       0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0       0      0      0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

Y 
16   000000000000000 
0000000000000000 

XX«-IFFT  Y 

XX 

1111111111111111 

0000000000000000 
XX+-XXC0;] 
XX 

1111111111111111 
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vpsEtnnv 

CO]   INW«-ZP  PSE  DATA;aiO;N:XN;XW 

CI]   A 

[2  3   A   COMPUTING  THE  POWER  SPECTRUM  ESTIMATION  CPERIODOGRAM3 

[33   A 

[4  3    QIO<-0 

[5  3    N*-CpDATA3+ZP 

[63    N«-2*CC2»N3 

[7  3    XN<-DATA.CCN-CpDATA3)pO) 

[8  3    XW^FFT  XN 

[93    INW<-C+i<CXW*2))*N 

[103   -»0 

[113  A 

[12  3  A   COMPUTING  THE  POWER  SPECTRUM  ESTIMATION  CPERIODOGRAM3  USING  FFT 

[13  3  A   FUNCTION 

[143  A 

[153  A   Y.KATZIR  .I.A.F. .AUGUST  1987 

[163  A 

[17  3    A      DATA<-THE   DATA  VECTOR    C1-D3 

[18  3    A      ZP«-LENGTH   OF   ZERO   PADDING   FOR   SMOOTH   PIX3T 

[193    A       INW.-THE   OUTPUT-SXCcj3 

[203    A 


VPSEB[Q3V 
[03   INBW<-SEG  PSEB  DATA:DIO;M:MM:XW;XN:L:K 
[13   A 

[2  3   A   COMPUTING  THE  POWER  SPECTRUM  ESTIMATION  USING  BARTIZTT  WINDOW 
[33   A 
[4  3    0I0«-0 
[53    M«-UCpDATA3*SEG) 
[63    XN«-CM.SEG)pO 
[7  3    ZP«-256-SEG 
[83    MM«-ClipXN)+ZP 
[93    XW*-CM.MM3pO 
[10  3   L«-K*-0 

[113  IX)OP1:XN[K;t.SEG3<-DATA[CSEGxK3+1.SEG3 
[12  3   -»(M>K«-K+13/LOOPl 
[133  LOOP2:XW[L;3«-ZP  PSE  XN[L;3 
[14]   -»CM>L<-L+13/IOOP2 
[15  3   INBW4-C+/XWD-5-M 
[16  3   ->0 
[173  A 

[183  A   Y.KATZIR, I.A.F. .AUGUST  1987 
[193  A 

[20  3  A   SEG*-LENGTH  OF  SEGMENTS 
[213  A   DATA«-THE  DATA  VECTOR  CI-D) 
[22  3  A   USING  THE  PSE  FUNCTION  C+FFT). 
[23  3  A   THE  FUNCTION  PADS  THE  SEGMENT  WITH  ZEROS  UP  TO  256  POINTS. 


[243  A 
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5CCCCCCCC 


c     c     c     c 

J    J    J    J 


5  3 


f«-240  PSZ  X 

Z*128tZ 

Y«al29 


Plot     of     Z    vs     Y 


25 
20 
15 

lO 
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■  ,,r  ■  ■ 
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1     ■    . 

si 
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I 

iLi 
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. —     j 

\ 

SO 


SO  90 

Y 


120 


1  SO 


5  5  5 


256 


555555555555555555 

rccceeecccceccccer 

*    m   m    c 
3   3    3   3 


5  5  5 


c     c      =■      c 


c      e      e      r      c      c 


J      -j      J      J 


c  c  c  e  e  c 

J  %f  Ml  */  >J  J 

c  c  c  r  =■  e 

>J  J  J  W  J  J 

c  c  c  c  c  e 

j  ^  j  -■  -<  .■ 


c     e     c     c     c     e- 


Z<-32  ?3ZB  X 
<>Z 

Zfl28tZ 


iOO 

60 
40 
20 


P 1  o-t     of     2    s^»     V 

<  4     se-srmen-t     of     32     points.  ) 


l   ■     ■     i     •     *    r  - « — « — r    i — ■ — i — ■ — "— : 
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n 
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150 
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VXSQRTCD3V 

CO]  XZ«-XSQRT  XR 

CI]  B 

C2]  «  COMPLEX  SQUARE  ROOT  CXR*0.5) 

C3]  A 

C4]  X2<-C0.5-Oa.>XRDx(2x(C+/XRxXR)*0.5D  +  I  1  0  +.xXRD*0.5 

C5]  XZ«-COa.>  l  0  /XR)eXZ,CaiO-0.5](0  1  +.  xXR)-i-2xXZ+XZ  =  0 


X 

0.5  0.866 

XSQRT  X 
0.86602  0.49999 


VXTIMECQ]7 

CO]  XZ+-XL  XTIME  XR 

CI]  A 

C2]  A  COMPLEX  TIMES  CXLxXR) 

C3]  fl 

C4]  XZ«-CC1  0  /XL)x  1  0  /XR3-C0  1  /XL)x  0  1  /XR 

C5]  XZ«-XZ,CDI0](C1  0  /XLJx  0  1  /XRD  +  CO  1  /XL)x  1  0  /XR 


X 

0.5  0.866 

X  XTIME  X 
"0.49996  0\866 
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VXMAGNCD37 
CO]       Z«-XMAGN   XR 
CI]       A 

[23       fl    COMPLEX  MAGNITUDE    (IXR)    —RETURN    IS    REAL 
C3]       R 
C4]         Z«-C+r<XR*2)*0.5 


X 

0.5    0.866 

XMAGN   X 
0.99998 


VXPHASCD3V 

CO]  Z«-XPHAS    XR 

CI]  A 

C2]  fl    COMPLEX   PHASE   ANGLE   —RETURN   IS   REAL 

C3]  P 

C4]         Z«-   0    1    +.XXR 
C5]         Z<-(xZ  +  Z  =  0)x(0v.j«XR)x-2oCl    0    +.  xXR)-K0a.  =  XR)  +  O/XRxXRD*0.5 


X 

0.5    0.866 


XPHAS    X 


1.0472 
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VXC0NJCD3V 
C03  XZ«-XCONJ  XR 
CI]   R 

C2]       ft    COMPLEX   CONJUGATE    C+XR) 
C31       ft 
C4]         XZ<-C1    0    7<XFJ.CniO]-    0    1    /XR 


X 

0.5    0.866 

XCONJ   X 
0.5    "0.866 


7XEXP0CQ]? 
CO]  XZ<-XEXPO  XR 
CI]       ft 

C2]       ft    COMPLEX   EXPONENTIAL    (*XR) 
C3]       ft 
C4]         XZ«-CCpXR)p*    1    0    ><XR)x    2    1    o.o    0    1    -t-.xXR 


X 

0.5  0.866 

XEXPO  X 
1.0682  1.2559 
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?PICD]V 
CO]   Y«-PI 
C1J    Y«-oi 


PI 
3.1416 


?COSCD]V 
CO]      Y<-COS    X 
CI]         Y<-2ox 


COS    CPI-5-4) 
0.70711 


vsiNccnv 

CO]      Y<-SIN   X 
CI]         Y*-lox 


SIN    CPI-r2) 
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) LOAD  0  EC4440 
0  EC4440  SAVED  2/23/1988  14:26:16 


NAVAL  POSTGRADUATE  SCHOOL  MONTEREY.  CALIFORNIA 

ECE-4440  WORKSPACE 


VERSION  1.0 


THESIS  BY 
Y.KATZIR.I.A.F. 
ADVISOR: PROF.  C.W.  THERRIEN.  , 

SEPTEMBER  1987 


)FNS 

BARTLETTC        BARTLETTR        CC       CC2D     CHEB     COEFF    COS      FCC 
FCC2D    FFT      FFT2D    FLCV     HAMMINGC         HAMMINGR         HANNINGC 
HANNINGR         IFFT     IFFT2D   ILPFILT  INFORMATION      LCV      MCCLEL   MEAN 
PI       PROTFILT         RECTANR  SHIFT    SHIFT2D  SIN      TRANSFNC         UNSHIFT 
UNSHIFT2D        XCONJ    XEXPO    XMAGN    XPHAS    XSQRT    XTIME 
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7CC2DCD3V 
C03   Y<-X  CC2D  H:I1:I2;DI0;N1;N2:T1;T2 
£13   a 

C2  3   A    2-D  CIRCULAR  CONVOLUTION 
C33   A 

C43    •+CCA/(pX)=pH)A2  =  ppX)  /BEGIN 
[53    □«-' IMPROPER  ARGUMENTS' 
[63    -*0 

[7  3   BEGIN:  DIO«-0 
[83    Nl«-lTpX 
[93    N2«-lipX 
[103   Y«-CN2.N13pJ) 
[113   Il«-0 

[123  LOOP1  :T1<-C-C  11+13  34*>H 
£13  3   I2*-0 

[143  LOOP2:T2«-C-CI2  +  13  3eeTl 
C153   Y[I2:I13«-+/  +  /XxT2 
[163   -KN2>I2«-I2  +  13/LOOP2 
C17  3   -KN1>I1<-I1+13/L00P1 
[13  3   -»0 
[193  A 

C203  A   THE  INPUTS  NEED  TO  BE  2-D  AND  THE  SAME  SIZE. 
C213  A 

C223  A   Y.KATZIR  I.A.F. .SEPTEMBER  1987 
[233  A 
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X 

2  2  10 
110  0 
0  0  0  0 
0    0    0    0 

H 
110    0 
110    0 
10   0   0 
0    0    0    0 

Y«-X   CC2D  H 

Y 

2  4    3    1 

3  6  4  1 
3  4  2  0 
110    0 


Zl 
1   2    3    4    5    6    7    8    9    10    11    12    13    14    15    16 

X   CC2D   21 

IMPROPER  ARGUMENTS 


Z2 


1  2 
3  4 


X  CC2D  Z2 

IMPROPER  ARGUMENTS 
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7FFT2DCQ37 

[03  Y«-FFT2D  X;QIO;NROW:NCOL;I:J;G 

CI]  fl 

[2  3  A   2-DIMENSIONAL  FFT 

C33  A 

[4  3  QIO«-0 

[53  A   CONVERT  REAL  ARRAY  TO  COMPLEX  ARRAY 

[63  I«-pCpX) 

[73  -KI  =  3  3/CONTIN 

[8  3  NROW<-pX[;13 

[93  NCOL*-pX[l:J 

[103  G«-C2.NROW.NCOL3pO 

[113  G[0::3«-x 

[123   X«-G 

[13  3  A  COMPUTATION  OF  FFT 

[143  CONTIN:NROW«-pX[l:;13 

[15  3   NCOL«-pX[l:l;3 

[16  3   G*-XxO 

[173   Y«-G 

[  18  3   H-J<-0 

[193  COLLOOP:G[;:J3«-FFT  X[  :  ;J3 

[20  3   -»CNCOL>J«-J*13/COLLOOP 

[213  ROWLOOP:Y[;l;3«-FFT  G[;I;3 

[22  3   -KNROW>I«-I  +  l)/ROWLOOP 

[23  3  A   IF  OUTPUT  ARRAY  IS  REAL  VALUE. 2ND  PLANE  IS  DELETED. 

[243   G[0;;3<-Y[0;;3+Y[l;;3 

[253   G[l;;3«-Y[l::3 

[263   !*+/♦/+/  CG*YJ 

[273   -KI  =  0)/STPNOW 

[283   •♦0  » 

[293    STPNOW:Y<-C1,NROW,NCOLJTY 

[303      Y«-CNROW,NCOL)pY 

[313    A 

[32  3  A   FFT  FUNCTION  HAVE  BEEN  USED. 

[333  A 

[343  A   Y. KATZIR. I. A. F. .SEPTEMBER  1987 

[353  A 
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VIFFT2DCDJV 
C03   Y<-IFFT2D  X:D.IO:NROW;NCOL;I;J;G 
CI]   A 

C2]   A   2-DIMENSIONAL  INVERSE  FFT 
C33   B 
C4  3    010<-0 

C5  3   A   CONVERTING  REAL  ARRAY  TO  COMPLEX  ARRAY 
C6]    I<-pCpX3 
C7]    -»CI  =  3)/CONTIN 
C8]    NROW*-pXC;i: 
C9]    NCOL«-pXCl;3 
C103   G«-C2.NROW.NCOL3p0 
till   GCO:;]«-X 
C12]   X«-G 

C13]  A  COMPUTATION  OF  IFFT 
C14]  CONTIN:NROW«-pXC0:;l] 
C15]   NCOL«-pXC0;l;] 
C16]   G«-XxO 
C17]   Y«-G 
C18]   I<-J«-0 

C19]  COLLOOP:GC;:J]*-IFFT  XC  : :  J] 
C20]   -»CNCOL>J«-J+13/COLLOOP 
C21]  ROWLOOP:YC;I;]<-IFFT  GC;I:] 
C22]   -*CNROW>I<-I+l)/ROWLOOP 

C23]  R   IF  OUTPUT  ARRAY  IS  REAL  VALUE. 2NT  PLANE  IS  DELETED. 
C241  GC0;:]*-YC0;;]+YC1;;] 
C25]   GC1;;]<-YC1;;] 
C26]   I*-+/*/  +  /CG9<Y) 
C27]   -KI  =  0)/STPNOW 
C28]   -»0  , 

C29]  STPNOW:Y«-Cl„NROW.NCOL]TY 
C30]   Y*CNROW,NCOL3pY 
C31]  A 

C32]  a   IFFT  FUNCTION  HAS  BEEN  USED. 
C33]  A 

C34]  A   Y. XATZIR. I. A. F. .SEPTEMBER  1987 
C35]  A 


110 


HI 
10000000 
00000000 

oooooooo 

00000000 

oooooooo 
oooooooo 
oooooooo 
oooooooo 

Y«-FFT2D  HI 


11111111 
11111111 
11111111 
11111111 
11111111 
11111111 
11111111 
11111111 

H2«-IFFT2D  Y 

H2 
10000000 
OOOOOOOO 
OOOOOOOO 

oooooooo 
oooooooo 
oooooooo 
oooooooo 
oooooooo 

XI 

2    1 
1    0 

1  1 

1    0 


Y+-FFT2D   XI 
Y 


4      2 

2  0 

3  1 
1    "I 


X2*-IFFT2D   Y 
X2 


2  1 

1  0 

1  1 

1  0 


111 


1110000000000000 
1110000000000000 
1110000000000000 
0000000000000000 
0000000000000000 
0000000000000000 
0000000000000000 
0000000000000000 
0000000000000000 
0000000000000000 
0000000000000000 
0000000000000000 
0000000000000000 
0000000000000000 
0000000000000000 
0000000000000000 

ZK-FFT2D  Z 
Z2*-XMAGN  Zl 

Z2 

9      8.5  7.2  5.3  3  0.7    1.2    2.5    3      2.5    1.2 

0.7  3  5.3  7.2  8.5 

8.5    8.1  6.9  5  2.8  0.67   1.2    2.4    2.8    2.4    1.2 

0.67  2.8  5  6.9  8.1 

7.2  6.9  5.8  4.3  2.4  0.57   1      2      2.4    2      1 
0.57  2.4  4.3  5.8  6.9 

5.3  5  4.3  3.1  1.8  0.41   0.73   1.5    1.8    1.5    0.73 
0.41  1.8  3.1  4.3  5 

3      2.8  2.4  1.8  1  0.23   0.41   0.85   1      0.85   0.41 

0.23  1  1.8  2.4  2.8 

0.7    0.67  0.57  0.41  0.23  0.055  0.097  0.2    0.23   0.2    0.097 

0.055  0.23  0.41  0.57  0.67 

1.2    1.2  1  0.73  0.41  0.097  0.17   0.35   0.41   0.35   0.17 

0.097  0.41  0.73  1  1.2 

2.5    2.4  2  1.5  0.85  0.2    0.35   0.72   0.85   0.72   0.35 

0.2  0.85  1.5  2  2.4 

3      2.8  2.4  1.8  1  0.23   0.41   0.85   1      0.85   0.41 

0.23  1  1.8  2.4  2.8 

2.5    2.4  2  1.5  0.85  0.2    0.35   0.72   0.85   0.72   0.35 

0.2,  0.85  1.5  2  2.4 

1.2  1.2  1  0.73  0.41  0.097  0.17   0.35   0.41   0.35   0.17 
0.097  0.41  0.73  1  1.2 

0.7    0.67  0.57  0.41  0.23  0.055  0.097  0.2    0.23   0.2    0.097 

0.055  0.23  0.41  0.57  0.67 

3      2.8  2.4  1.8  1  0.23   0.41   0.85   1      0.85   0.41 

0.23  1  1.8  2.4  2.8 

5.3  5  4.3  3.1  1.8  0.41   0.73   1.5    1.8    1.5    0.73 
0.41  1.8  3.1  4.3  5 

7.2    6.9  5.8  4.3  2.4  0.57   1      2      2.4    2      1 

0.57  2.4  4.3  5.8  6.9 

8.5    8.1  6.9  5  2.8  0.67   1.2    2.4    2.8    2.4    1.2 

0.67  2.8  5  6.9  8.1 
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i^IillS^i-dc^o) 


output  rue.,/ 


I!!i!£^y domino 


OUTPUT  FNC.U regency  domain) 
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7SHIFTCI33V 

CO]  Y«-SHIFT  X 

CI]  A 

C23  A   SHIFTING  PERIODIC  DATA  TO  AN  INTERVAL  STARTING  AT  THE  ORIGIN 

[33  A 

C43    Y«-CLC~l+pX3*2)<t>X 
C53    -»0 

[63  A 

C7]  A   SOURCE :M-D  SIGNAL  PROCESSING  W.S. 

C8]  A 


X 

000010000 


X1«-SHIFT  X 

XI 
1000  0-  0000 


VUNSHIFTCD3V 

C03  Y<-UNSHIFT  X 

C13  A 

C23  A   SHIFTING  PERIODIC  DATA  TO  AN  INTERVAL  CENTERED  AT  THE  ORIGIN 

C33  A 

C4]  Y<-C-LC-l+pX)*2)<J>X 

C53  -»0 

C63  A 

C73  A      SOURCE:M-D   SIGNAL   PROCESSING   W.S. 

C83  A 


XI 

100000000 

X2<-UNSHIFT   XI 

X2 
000010000 
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VSHIFT2DCD3V 

CO]  Y<-SHIFT2D  X 

CI]  ft 

[23  A   SHIFTING  2-D  PERIODIC  DATA  TO  AN  INTERVAL  STARTING  AT  THE  ORIGIN 

[3]  A 

C4  3  Y«-CLC-H-lTUpX)-r2)<J>ClC-l+lfpX3-i-2)eX 

[53  ->0 

£6]  A 

[73  A   SOURCE:M-D  SIGNAL  PROCESSING  W.S. 

[83  A 


7UNSHIFT2D[Q3V 

[03  Y«-UNSHIFT2D  X 

[13  A 

[2  3  A   SHIFTING  2-D  PERIODIC  DATA  TO  AN  INTERVAL  CENTERED  AT  THE  ORIGIN 

[33  A 

[4  3    Y<-C-LC-l+lTllpX)*23*C-lC-l+lTpX3*2)eX 
[53    -»0 

[S3  A 

[73  A   SOURCE: M-D  SIGNAL  PROCESSING  W.S. 

[83  A 
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000000000 
000000000 
000000000 
000000000 
000010000 
000000000 
000000000 
000000000 
000000000 

ZK-SHIFT2D  Z 

Zl 
100000000 
000000000 
000000000 
000000000 
000000000 
000000000 
000000000 
000000000 
000000000 

Z2«-UNSHIFT2D  Zl 

Z2 
000000000 
000000000 
000000000 
000000000 
000010000 
000000000 
000000000 
000000000 
000000000 
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7ILPFILTCD]? 

CO]  W«-R  ILPFILT  P;QIO:I:J;RR;RB;PB 

CI]  ft 

C2]  ft   IDEAL  LP  FILTER  C CIRCULAR  BASE)  FREQUENCY  RESPONSE 

C3]  ft 
C4]    DIO-0 

C5]  W<-CP.P)pO 
C  6  ]    I«-J«-0 

C7]  PB«-(P+2  3-l 

C8]  RB«-(PxR)  +  C2xPI) 

C9]  LOOP:RR«-CCCI-PB)*2]  +  (CJ-PB)*2))*0.5 

CIO]  -»(RR>RB3/ZERO 

Cll]  WCI;J]*-1 

C12]  ZERO  :-KP>J<-J+l) /LOOP 

C13]  J<-0 

C14]  -KP>I<-I  +  l)/LOOP 

[153  -»0 

[16]  A 

[17]  ft   R«-THE  RADIUS  OF  THE  BASE  CGAIN=13 

[18  3  ft   P«-THE  DIMENSION  OF  THE  MATRIX 

[19]  ft   W:THE  IDEAL  FREQUENCY  RESPONSE 

[20]  ft 

[21]  ft   Y.KATZIR.I.A.F. .  OCTOBER  1987 

[22]  ft 
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VHAMMINGRC03V 
CO]      WR<-II   HAMMINGR  N:OIO:HW1;HW2  :  J;I  ;ROW;COL;S 
CI]       ft 

C2  3       A      RECTANGULAR   BASE   HAMMING   WINDOW 
C33       fi 
C4  3         QIO«-0 
C53         WR<-(N.N)pO 
C63         HWl«-HW2«-IIpO 
C73         S«-CII-1)*2 
C  8  ]         ROW<-0 
C93      ROWLOOP:  COL<-0 
CIO  3       I<-ROW+C2xSD 

C113      HWlCROW]<-0.54  +  C0.46xCOSCCPIxCROW-S)D-s-S)) 
C12  3    COLLOOP:J<-COL+C2xS) 

C133      HW2CCOL3«-0.54  +  C0.46xCOSC  CPIx(COL-S3  3-s-SD) 
C14]      WRCI:J3<-HWlCROW]xHW2CCOL3 
C15]      -KII>COI>COL+1)/COLLOOP 
C16]       -KII>ROW«-ROW-i-l] /ROWLOOP 
C17]       -»0 
C18]     A 

[193    A       II+-WINDOW   BASE    DIMENSION 
C20  3    ft      N«-MATRIX   DIMENSION 
C213    A      WR:THE   WINDOW  MATRIX 
[22]    ft 

C233    A      Y. KATZIR. I. A. F. .SEPTEMBER    1987 
C243    A 
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N*32 

UI+I1    HAMMINGS  N 
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VHAMMINGCCCnv 
[03   WOII  HAMMINGC  N:DI0;I ; J;PR;RB;RR 
CI]   A 

[2  3   A   CIRCULAR  BASE  HAMMING  WINDOW 
[33   A 
[43    DIO<-0 
[53    WOCN.N)p0 
[  6  3    I<-J<-0 
[7  3    PR«-CN*2  3-1 
[83    RB«-CIIxN3  +  C2xPI3 

[93   I£OP:RR«-(CCI-PR)*2)  +  CCJ-PR3*2  3D*0.5 
[103   -»CRB<RR3/C0NT 

[113   WC[J:I3<-0.54  +  C0.46xCOSCCPIxRR)-rRB3) 
[123  C0NT:-»CN>J«-J+13/L00P 
[13  3   J*-0 

[143   -KN>I<-3>13/I£OP 
[153   -+0 
[163  A 

[17  3  A   IK-THE  RADIUS  OF  THE  BASE 
[18  3  A   N<-THE  DIMENSION  OF  THE  MATRIX 
[19  3  A   WC:THE  HAMMING  WINDOW 
[203  A 

[213  A   Y.KATZIR  .I.A.F. .  OCTOBER  1987 
[223  A 
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II«-0.4x?I 
N*32 

UOII    HAMMINGC  N 
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VRANNINGRCD3G 
CO]      WR«-H   HANNINGR  N:DIO;HWl;HW2 ; J;I;ROW;COL;S 
C13       A 

C2  3       A      RECTANGULAR   BASE   BANNING   WINDOW 
C33      a 
C4  3         DIO<-0 
CS3         WR«-CN.N3pO 
C63         HWl«-HW2«-IIpO 
C73         S*CII-13*2 
C  8  3         ROW<-0 
C9  3      ROWLOOP:CO3>0 
CIO  3       I«-ROW+C2xS3 

C113      HW1CROW3-K3. 5xCl+CCOSC  CPIxCROW-S3  3+S3  3  3 
C12  3    COUJOOP:J*-COL*C2xS3 

C13  3      HW2CCOL3*-0.5x(l+(COSCCPIx(COL-SD)H-S3  3  3 
C143      WRCI:J3«-HWlCROW3xHW2CCOL3 
C153      -KII>COL«-C03>13/C01iOOP 
C163      -KII>ROW*-ROW+13/ROWI£OP 
C173      -+0 
C183    A 

C19  3    A      IK-WINDOW   BASE   DIMENSION 
C20  3    A      N<-MATRIX   DIMENSION 
C213    A      WR:THE  WINDOW  MATRIX 
C223    A 

C233    A      Y. KATZIR, I. A. F. .SEPTEMBER   1987 
C243    A 
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VHANNINGCED3V 
CO]      WOII   HANNINGC   H;DXO;I : J:PR;RB;RR 
£13       A 

£2  3       A      CIRCULAR   BASE   HANNING  WINDOW 
£33      n 
£4  3         0IO<-0 
£53         WC*-CN.N3pO 
£  6  3         I<-J<-0 
£7  3         PR<-CN+2)-l 
£83         RB<-CIIxN3-K2xPI] 

£93      I£OP:RR<-CCCI-PR3*2)  +  CCJ-PR)*2]  )*0.5 
£103      -KRB<RR)/CONT 

£113      WCEJ:I3<-0.5xCl+CCOSCCPIxRR)-rR£O)) 
£123    C0NT:-KN>J*-J+1]/I£0P 
£13  3      J<-0 

£143   -»CN>I<-I+1)/L00P 
£153   -»0 
£163  n 

£17  3  A   II<-THE  RADIUS  OF  THE  BASE 
£18  3  A   N<-THE  DIMENSION  OF  THE  MATRIX 
£19]  A   WC:THE  HANNING  WINDOW 
£203  a 

£213  A   Y.KATZIR.I.A.F.  .  OCTOBER  1987 
£223  A 
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VBARTLETTR C 0 3 7 
CO]   WR<-II  BARTLETTR  N:DI0;HW1;HW2  ;  J;I  :ROW;COL;S 
C1J   ft 

C2J   A   RECTANGULAR  BASE  BARTLETT  WINDOW 
C3]   A 
C4]    QIO<-0 
C53    WR<-(N.NDpO 
[6]    HWl«-HW2*-IIpO 
C7]    S«-CII-l)+2 
C  8  3    ROW4-0 
C9]   ROWLOOP:C0L«-0 
C103   I<-ROW+(2xS) 
C113   HW1  CROW] «-l- CI  CROW-SDKS 
C12]  COLLOOP:J*-COL+C2xS) 
C13J   HW2CC0U<-1-C  I  CC0L-S3DH-S 
C143   WRCI;J]«-HWlCROW3xHW2CCOLJ 
C15]   -KII>COL<-COL+l)/COLLOOP 
[16  3   -KII>ROW«-ROW+lD/ROWLOOP 
C173  +0 
C18]  A 

C19  3  A   II«-WINDOW  BASE  DIMENSION 
C20  3  A   N«-MATRIX  DIMENSION 
C213  A   WR:THE  WINDOW  MATRIX 
C223  A 

C233  A   Y. KATZIR. I. A. F. .SEPTEMBER  1987 
C243  A 
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VBARTLETTCCD3V 
CO]      WOII    BARTLETTC  N:DIO;I:J;PR:RB;RR 
[1]       R 

[23       «      CIRCULAR   BASE   BARTIETT  WINDOW 
[3]       ft 
C4]         DIOK) 
[53         WC«-CN.NDpO 
C63         I«-J<-0 
[73         PR«-CN*2D-1 
[83         RB«-CIIxN)-K2xPI) 

£9  3       IOOP:RR«-CCCI-PR)*2)  +  CCJ-PR)*2n*0.5 
C103      -KRB<RRD/CONT 
[113      WC[J:I3<-1-CRR*RB) 
[12  3    C0NT:-KN>J«-J+13/L00P 
[13  3      J<-0 

[143   -»(N>3>I  +  1)/L00P 
[153   -»0 
[163  ft 

[17  3  ft   II«-THE  RADIOS  OF  THE  BASE 
[18  3  A   N«-THE  DIMENSION  OF  THE  MATRIX 
[19  3  fl   WC:THE  HAMMING  WINDOW 
[203  fl 

[213  B   Y.KATZIR  .I.A.F. ,  OCTOBER  1987 
[223  a 
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VRECTANRCCHV 
CO]      WR^II   RECTANR  N:DIO:HWl:HW2 : J;I ;ROW;COL;S 
CI]       A 

C23       A      RECTANGULAR   BASE   RECTANGULAR   WINDOW 
C3]       A 
C4]         DIO«-0 
C5]         WR«-CN.N]p0 
C6]         HWl<-HW2<-IIpO 
C7]         S«-CII-1]*2 
C  8  ]         ROW«-0 
C9]      ROWLOOP:COL«-0 
CIO]      I<-ROW+C2xS) 
Cll]      HW1CR0W]<-1 
C12]    COLLOOP:J«-COL+C2xSD 
C13]      HW2CC0L]<-1 
C14]      WRCI:J]«-HWlCROW]xHW2CC0L] 
C15]      -»CII>C0L<-C0L+1]/C0LL00P 
C16]       -KII>ROW<-ROW+l]/ROWLOOP 
C17]      -»0 
[131    A 

C19]    A      II«-WINDOW   BASE    DIMENSION 
C20]    fl      N<-MATRIX   DIMENSION 
C21]    A      WR:THE   WINDOW  MATRIX 
C22]    A 

C23]    A      Y. KATZIR. I. A. F. .SEPTEMBER   1987 
C24]    A 
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CO]  HPROT«-PROTFILT  SIZE : 010 ; W :N1 ; N2 :N 

C1J  A 

C2]  A   PRODUCING  A  LOWPASS  PROTOTYPE  FILTER 

C3]  A 

C4]  oio<-m«-o 

C51    N1*-SIZE*2 

C6]    N2«-N1-1 

C7J    HPROT«-SIZEoO 

C8]  L00P:W«-CCN-N2DxPI)-rNl 

C9]    W4-IW 

C10I   -»(WX0.4xPI))/ZERO 

Cll]   HPR0TCN]«-1 

C12]  ZERO:  •*(SIZE>N«-N*1) /LOOP 

C13J   HPROT«-SHIFT  HPROT 

C14J   -*0 

CI53  A 

C16J  A   SIZE«-THE  SIZE  OF  THE  PROTOTYPE  FILTER  (NO.  OF  SAMPLES) 

C17]  A   THE  SHIFT  FUNCTION  HAS  BEEN  USED. 

CIS]  A 

C19]  A   Y.KATZIR.I.A.F. .  OCTOBER  1987 

C20I  A 
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VCOEFFCC33? 

CO]  a«-COEFF  Hp;aiO:SI2E:HH:fl 

CI]  a 

C2  3  »   COMPUTING  THE  COEFFICIENTS  FROM  THE  PROTOTYPE  FILTER  (USING  IFFT) 

C3]  A 

C4]  aio«-o 

C5]  SIZE«-(pHp)+2 

C6]  a«-SIZEp0 

C7]  H«-IFFT  Hp 

C8]  HHMI.SIZEJTH 

C9]  HH«-SIZE<»HH 

CIO]  aClSIZE]«-2xHHCl.SIZE] 

Cll]  aC0]«-aC0]+2 

C12]  -»0 

C13]  a 

C14]  ft   Hp«-THE  PROTOTYPE  SAMPLES 

C15]  a   a: THE  IMPULSE  RESPONSE  COEFFICIENTS 

C16]  «   THE  IFFT  FUNCTION  HAS  BEEN  USED. 

[17]  » 

C13]  »   tf.KATZIR.I.A.F.  .  OCTOBER  1987 

C19]  fl 
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7TRANSFNCCD3V 

C03  Fw«-TRANSFNC  SIZE:C]I0:A:B;R0W:C0L:W1:W2:N1:N2 

C13  a 

C23  A   PRODUCING  2-D  TRANSFORMATION  FUNCTION  C EXAMPLE) 

C33  A 

C43  0IO«-0 

C53  Fw«-<  SIZE,  SIZE  )p0 

C6]  ROW«-0 

C7J  Nl«-SIZE+2 

C83  N2«-N1-1 

C9]  ROWLOOP:Wl<-(CROW-N2)xPI)+HI 

CIO 3  A«-COSCWl) 

C113  COL+-0 

C 12  3  COLLOOP: W2*< ( COL-N2 ) *PI ) +N1 

C13  3  B«-COSCW2) 

C143  FwCROW:COL3«-(C-l)«-A+B+CAxB))+2 

CIS  3  -»(SIZE>COL«-COL*l)  /COLLOOP 

CIS  3  -»CSIZE>ROW«-ROW*l)/ROWLOOP 

C17  3  -»0 

C181  A 

C193  A   THE  FUNCTION  PRODUCING  THE  TRANSFORMATION  MATRIX  IN  ORDER  TO 

C203  A   DESIGN  FIR  FILTER  USING  McCLELLAN  TRANSFORMATION 

C213  A   SIZE«-NO.  OF  SAMPLES  IN  ul.  C<u2  WILL  BE  THE  SAME  SIZE  AUTOMATICALLY) 

C22  3  A   FvwtTHE  TRANSFORMATION  MATRIX 

C233  A   IN  ORDER  TO  CHANGE   FC«1.«2) .REPLACE  LINES  C103.C133.  AND  C14) 

C24)  A   WITH  THE  NEW  CHOICES. 

C253  A 

C263  A   Y.KATZIR.I.A.F. .  OCTOBER  1987 

C273  A 
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VTRANSFNCtn:* 
C03   Fvw«-TRANSFNC  SIZE:aiO:A:B:ROW;COL:Wl;W2:Nl:N2 

ci:     a 

C23       A      PRODUCING   2-D   TRANSFORMATION   FUNCTION    (EXAMPLE) 

C33       A 

C4  3         0IO«-0 

C53         Fw«-CSIZE.SIZE)pO 

C63         ROW4-0 

C73         NK-SIZE+2 

t83         N2«-N1-1 

C93      ROWLOOP:Wl«-CCROW-N23*PI)+Nl 

C103      A«-SINCW1) 

C113      COI>0 

C123    COLLOOP:W2«-(CCOL-N2)xPI)+Nl 

C13  3      B«-SINCW2) 

C143      FwCROW;COL3«-AxB 

C15  3   -+(SIZE>COL«-COL*l)/COLLOOP 

C163   -»CSIZE>ROW«-ROW*l)/ROWLOOP 

[173   -»0 

[133  A 

[19  3  R   THE  FUNCTION  PRODUCING  THE  TRANSFORMATION  MATRIX  IN  ORDER  TO 

C203  A   DESIGN  FIR  FILTER  USING  McCLELLAN  TRANSFORMATION 

[213  A   SIZE«-NO.  OF  SAMPLES  IN  ul.  Cu2  WILL  BE  THE  SAME  SIZE  AUTOMATICALLY) 

C223  A   Fw:THE  TRANSFORMATION  MATRIX 

[233  A   IN  ORDER  TO  CHANGE   F C ul , u2 J . REPLACE  LINES  [103. [133.  AND  [143 

C243  A   WITH  THE  NEW  CHOICES. 

C253  A 

C2S3    A      Y.KATZIR.I.A.F.  .    OCTOBER   1987 

[273    A 
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7CHEBCD3V 

CO]  Y«-N  CHEB  X 

C13  fl 

C2  3  A   EVALUATING  CHEBYSHEV  POLYNOMIALS 

C3J  A 

[43    ->CCN0).N*1)/A.B 
[53    Y«-C2xXxCN-l)CHEB  X3-CN-23CHEB  X 
C63    -»0 

C73  A:Y«-CpX)pl 
C33    ->0 

C93  B:Y«-X 
[103   -»0 

£1X3  A 

[12  3  A   N«-THE  POLYNOMIALS  ORDER 

[13  3  A   X«-THE  FUNCTION  TO  BE  EVALUATED 

[143  A 

C153  A   SOURCE:MDSP  M.S. 

C163  A 


X 

12  3 

0  CHEB  X 
111 

1  CHEB  X 
12  3 

2  CHEB  X 
1  7  17 

3  CHEB  X 
1  26  99 

4  CHEB  X 
1  97  577 

5  CHEB  X 
1  362  3363 

6  CHEB  X 
1  1351  19601 
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^MCCLELCtm 

CO]  Hwv«-a   MCCLEL  F;aiO:N:INTER;NUMBER:INTERl:INTERO 

CI]  fl 

C2]  A      DESIGNING    2-0    FILTER  USING   MEcCLELLAN   TRANSFORMATION 

C3]  A 

C4]      aio«-o 

C5]    INTERO«-INTERl«-INTER«-Hww«-<pF]pO 

C6]    NUMBER«-pa 

C7]    N«-0 

C8]  LOOP:-KCN=0).N=1)/A.B 

C9]    INTER«-C2xFxINTERl)-INTER0 

CIO]  -*C 

Cll]  A:INTER«-CpF)pl 

C 12  ]   INTER1*-  INTER 

ci3]  -k: 

C14  1  B:INTER«-F 

CIS]  C:a«-' CALCULATED  CHEB.  POLYNOM  OF  ORDER: '.IN 
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APPENDIX  D 
COMPUTER  ASSIGNMENTS  AND  SOLUTIONS 

This  appendix  provides  a  few  samples  of  computer  assignments  and  related 
solutions  using  the  software  package.  These  are  in  addition  to  the  examples  in 
Chapter  III.  The  computer  assignments  in  EC  3400  were  taken  from  FIRST  PRIN- 
CIPLES OF  DISCRETE  SYSTEMS  AND  DIGITAL  SIGNAL  PROCESSING  by 
Professors  R.  D.  Strum  and  D.  E.  Kirk  [4].  The  computer  assignments  in  EC  3410 
and  EC  4440  were  taken  from  Professor  C.  W.  Therrien's  homework  assignments. 


EC  3400  COMPUTER  ASSIGNMENT  [4] 

The  purpose  of  a  bandstop  filter  is  to  remove  the  effects  of  the  vibrations  in 
the  frequency  range  6  Hz  to  15  Hz  from  a  control  loop.  The  sampling  frequency  is 
fs=100  Hz. 

a.  Using  the  Fourier  Series  design  procedure,  calculate  the  filter  coefficients  with 
a  rectangular  window. 

b.  Plot  the  magnitude  portion  of  the  frequency  response  for  a  filter  having  thirty- 
one  coefficients. 

c.  Repeat  part  (b)  using  a  von-Hann  window. 


SOLUTION 

1.  By  using  the  function  DIGFREQ  the  digital  upper  and  lower  frequencies  are 
calculated  (Figure  D.l).  Figure  D.2  describes  the  ideal  frequency  response  to 
be  generated. 
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Figure  D.l.    Calculating  the  Digital  Frequency  Using  the  Function  DIGFREQ. 


H 


1 

o.  a 
o.  6 

0.4 

O.  2    f- 
O 


Li 


20 


C9J 


40  60 

Cx?I/64] 


SO 


Figure  D.2.    The  Ideal  Frequency  Response. 

2.  By  using  the  function  BSCOEFF  31  causal  impulse  response  coefficients  are 
calculated  (See  Figures  D.3  and  D.4). 

3.  The  function  HANNING  generates  31  von-Hann  window  samples  that  are 
multiplied  by  the  filter  coefficients.  The  result  is  windowed  coefficients  (See 
Figures  D.5  and  D.6). 

4.  By  using  the  function  FREQRES  the  frequency  responses  with  and  without 
window  are  generated.  The  following  shows  the  APL  commands  which  gener- 
ate the  frequency  responses,  and  Figures  D.7  and  D.8  show  the  plotted  results. 

HBS<-128  FREQRES  hBS 
HBSW<-128  FREQRES  hBSW 


143 


hBS«-31   BSCOEFF  TUL 

hBS 
'0.033694    "0.032561    "0.016485    "3.2342E~4    "1.0213E~3    "0.01371    "0.037409 

"0.032354    7.3543E~3    0.07206    0.12421    0.12619    0.063217    "0.042416    "0.14 
0.82    "0.14034    "0.042416    0.063217    0.12619    0.12421    0.07206    7.3543E"3 
"0.03285"4    "0.037409    "0.01371    "1.0218E~3    "3.2342E"4    "0.016435    "0.03 
"0.033694 


Figure  D.3.    Calculating  the  Impulse  Response  Coefficients. 
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Figure  D.4.    The  Impulse  Response  Coefficients. 
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Figure  D.5.    Generating  the  Hanning  Window. 
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Figure  D.6.    The  Hanning  Window. 
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Figure  D.7.    The  Frequency  Response  Magnitude. 
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Figure  D.8.    The  Frequency  Response  Magnitude  with  Window. 
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EC  3400  COMPUTER  ASSIGNMENT  [4] 

A  nonrecursive  bandpass  filter  is  to  be  designed  using  the  ID  FT  approach. 
The  ideal  characteristic  is  given  below: 
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Use  255  frequency  samples  and  find  51  filter  coefficients  using: 

a.  A  rectangular  window 

b.  A  Hamming  window 

c.  Plot  the  frequency  response  for  the  filters  of  parts  (a)  and  (b). 


SOLUTION 

1.  By  using  the  function  IDEALF  and  choosing  256  samples,  the  ideal  bandpass 
frequency  response  is  generated  (See  Figure  D.9). 

2.  By  using  the  function  FCOEFF  51  coefficients  are  calculated  (See  Figure 
D.10). 

3.  The  function  HAMMING  generates  51  Hamming  window  samples  that  are 
multiplied  by  the  filter  coefficients.  The  result  is  the  windowed  coefficients. 
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Figure  D.9.    Generating  the  Ideal  Frequency  Response  Using  255  Sam- 
ples. 
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Figure  D.10.    Calculating  the  Impulse  Response  Coefficients. 
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Figure  D.ll.    Generating  the  Hamming  Window  and  Calculating  the 
Frequency  Response. 
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Figure  D.12.    The  Frequency  Response  Plots. 
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Figure  D.13.    The  Frequency  Response  in  dB 
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Figure  D.14.    The  Frequency  Response  Phase. 
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EC  3410  COMPUTER  ASSIGNMENT  [5] 

On  a  floppy  disk  you  will  find  4  data  sets  called 

SLO.DAT 

SL1.DAT 

SL2.DAT 

SL3.DAT 
Each  of  these  represents  a  random  signal  with  512  time  points.  The  format  of  each 
file  is  identical.  You  can  look  at  the  file  with  an  editor.  The  first  line  contains  the 
number  512  displaying  the  number  of  points  in  the  data  set.  On  the  remaining 
lines,  having  four  numbers  per  line,  you  will  find  the  floating  point  values  of  the 
signal. 

1.  Transfer  these  files  to  your  APL  WS  using  the  UTILITY  function  GETDATA. 

2.  Plot  each  of  the  signals  versus  time  and  submit  the  plots  with  your  assign- 
ment. Tell  whether  you  think  each  of  these  signals  seems  to  be  more-or-less 
uncorrelated,  positively  correlated,  or  negatively  correlated. 

3.  Compute  the  mean  of  each  signal  by  averaging  in  time  and  write  down  your 
result  for  each  signal.  It  is  reasonable  based  on  your  plots? 

4.  Subtract  the  mean  from  each  signal  and  compute  and  plot  the  sample  corre- 
lation function.  Attach  the  plots  to  your  assignment.  Was  your  guess  about 
the  correlation  of  these  signals  in  part  2  correct? 


SOLUTION 

1.  Using  the  function  GETDATA  and  the  following  statements,  the  data  sets 
SLO,  SL1,  and  SL3  are  transferred  to  the  APL  workspace: 

SLO<-GETDATA  'SLO.DAT' 
SL1<-GETDATA  'SL1.DAT 
SL2 ^-GETDATA  'SL2.DAT 
SL3*-GETDATA    'SL3.DAT 
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2.  By  using  the  X-Y  plotting  procedure  in  STATGRAPHICS  the  data  plots  are 
generated.  (See  Figures  D.15,  D.16,  D.17,  and  D.18.)  Another  kind  of  plot 
can  be  obtained  by  using  the  STATGRAPHICS  Time  Series  Procedure.  The 
results  are  shown  in  Figures  D.19,  D.20,  D.21,  and  D.22. 
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Figure  D.15.    SLO  Data  Set. 
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Figure  D.16.    SL1  Data  Set. 
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Figure  D.17.    SL2  Data  Set. 
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Figure  D.18.    SL3  Data  Set. 

3.  From  the  plotting  of  each  data  set,  it  can  be  seen  that: 

a.  SLO  changes  randomly  with  little  consistency,  this  is  characteristic  of  a 
white  noise  sequence. 

b.  SL1  varies  more  slowly  and  stays  stable  for  a  while  with  positive  values 
and  then  changes  more-or-less  smoothly  to  negative  values.  This  indicates 
that  it  is  a  positively  correlated  data  set. 

c.  SL2  is  very  nearly  like  a  white  noise  sequence. 
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Figure  D.19.    SLO  Data  Set  Using  Time  Series  Procedure. 
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Figure  D.20.    SLl  Data  Set  Using  Time  Series  Procedure. 

d.  SL3  changes  at  almost  every  step  from  negative  to  positive.  This  is  char- 
acteristic of  a  negatively  correlated  data  set. 

4.  Calculate  the  mean  of  each  data  set  using  the  MEAN  function.  Then  subtract 
the  mean  of  each  data  set  from  the  data  set  in  order  to  achieve  a  zero  mean  data 
set.  Using  the  following  instructions,  calculate  the  first  150  autocorrelation 
lags  for  each  data  set: 

Rx:cl«-150  SAC?  SLO 
Rxx2*-150  SACF  SLl 
Rxx3*-150  SACF  SL2 
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Figure  D.21.    SL2  Data  Set  Using  Time  Series  Procedure. 
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Figure  D.22.    SL3  Data  Set  Using  Time  Series  Procedure. 

focx4<-150   SACF  SL3 

Using  the  X-Y  plotting  procedure  in  STATGRAPHICS  the  estimated  auto- 
correlation is  plotted  as  shown  in  Figures  D.23  to  D.26. 

Another  method  for  producing  the  autocorrelation  values  and  plot  them  is  to 
use  the  autocorrelation  procedure  in  the  TIME  SERIES  ANALYSIS  section  from 
the  STATGRAPHICS  Menu.  The  results  are  shown  in  Figures  D.27  to  D.30. 
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Figure  D.23.    Sample  Autocorrelation  Function  of  White  Noise. 
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Figure  D.24.    Sample  Autocorrelation  Function  of  Positive  Correlated 
Set. 

Looking  at  the  plots  verifies  that  our  characterization  of  the  data  sets  in  para- 
graph 3  was  correct. 
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Figure  D.25.    Sample  Autocorrelation  Function  of  Nearly  White  Noise. 
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Figure  D.26.    Sample  Autocorrelation  Function  of  Negative  Correlated 
Data  Set. 
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Figure  D.27.    Estimated  Autocorrelation  Using  Time  Series  Procedure. 
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Figure  D.28.    Estimated  Autocorrelation  Using  Time  Series  Procedure. 
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Figure  D.29.    Estimated  Autocorrelation  Using  Time  Series  Procedure. 
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Figure  D.30.    Estimated  Autocorrelation  Using  Time  Series  Procedure. 
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EC  3410  COMPUTER  ASSIGNMENT  [5] 

Compute  power  spectrum  estimates  for  the  data  sets  SLO,  SL1,  SL2,  and  SL3 
using  a  rectangular  window.  In  all  cases  zero-pad  your  data  appropriately  before 
giving  it  to  the  FFT  routine  so  that  you  will  have  computed  enough  points  of  the 
spectrum  to  give  a  smooth  plot.  (Plotting  say  256  points  of  the  spectrum  should 
be  sufficient).  Assume  the  data  has  been  sampled  at  a  spacing  of  T  =  0.1  ms. 
Label  the  frequency  axis  in  hertz.  Since  the  signals  are  real,  it  is  only  necessary  to 
plot  the  positive  frequencies. 

a.  Use  only  one  16-point  segment  (the  first  16  points)  of  each  data  set. 

b.  Use  32  16-point  segments. 

c.  Use  16  32-point  segments. 

d.  Use  4  128-point  segments. 
What  do  you  make  of  the  results? 


SOLUTION 

1.  The  first  step  is  to  extract  the  first  16  samples  of  each  data  set  and  compute 
the  power  spectrum  estimate  using  the  function  PSE.  In  order  to  obtain  a 
smooth  plot  the  data  sets  are  zero  padded  with  240  points.  (See  Figure  D.31.) 
The  power  spectrum  estimate  will  have  250  points  as  shown  on  the  plots  of 
Figures  D. 32  to  D. 35. 

2.  By  using  the  function  PSEB  (power  spectrum  estimation  with  Bartlett  win- 
dow), we  estimate  the  periodogram  of  the  SLO  data  set  using  30  segments  of 
17  points,  16  segments  of  32  points,  and  4  segments  of  128  points  as  seen  in 
Figure  D.36.  Figure  D.37  is  the  periodogram  plots  for  each  case.  The  above 
process  is  repeated  with  the  data  sets  SL1,  SL2,  and  SL3.  See  Figures  D.38 
to  D.43. 

3.  From  these  results  we  conclude  that  the  first  data  set  has  a  spectrum  like 
white  noise.    The  second  and  third  spectra  are  lowpass  with  bandwidths  of 
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Figure  D.31.    Calculating  the  Power  Spectrum  Estimates   Using  the 
Function  PSE. 

approximately  1  kHz  and  1.2  kHz.  The  fourth  spectrum  is  highpass  with 
bandwidth  of  approximately  1  Khz  Hz.  It  seems  that  the  best  estimate  is 
obtained  by  using  16  segments  of  32  points  each. 
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Figure  D.32.    The  PSE  of  SLO  Data  Set  Using  16  Points. 
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Figure  D.33.    The  PSE  of  SL1  Data  Set  Using  16  Points. 
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Figure  D.34.    The  PSE  of  SL2  Data  Set  Using  16  Points. 
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Figure  D.35.    The  PSE  of  SL3  Data  Set  Using  16  Points. 
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SL01«-17    PSEB   SLO 
SL01+-128TSL01 


SL02«-32    PSEB   SLO 
SL02«-12JiTSL02 


SL03«-128    PSEB   SLO 
SL03«-128TSL03 


Figure  D.36.    Calculating  Periodograms  of  SLO  Using  Bartlett  Win- 
dows. 
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Figure  D.37.    Periodograms  of  Different  Windows. 
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SL11«-17    PSEB   SL1 
SL1K-128TSL11 


SL12«-32    PSEB   SL1 
SL12«-128TSL12 


SL13«-128    PSEB   SL1 
SL13«-128TSL13 


Figure  D.38.    Calculating  Periodograms  of  SL1   Using  Bartlett  Win- 
dows. 
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Figure  D.39.    Periodograms  of  Different  Windows. 
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SL2 1<-17    PSEB   SL2 
SL21«-128TSL21 


SL22«-32    PSEB   SL2 
SL22«-128TSL22 


SL23«-128    PSEB   SL2 
SL23*-128TSL23 


Figure  D.40.    Calculating  Periodograms  of  SL2  Using  Bartlett  Win- 
dows. 
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Figure  D.41.    Periodograms  of  Different  Windows, 
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SL3 1<-17    PSEB   SL3 
SL31+-128TSL31 


SL32«-32    PSEB   SL3 
SL32«-128TSL32 


SL33<-128    PSEB   SL3 
SL33«-128TSL33 


Figure  D.42.    Calculating  Periodograms   of  SL3  Using  Bartlett   Win- 
dows. 
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Figure  D.43.    Periodograms  of  Different  Windows. 
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EC  4440  COMPUTER  ASSIGNMENT-FILTER  DESIGN  [6] 

This  problem  deals  with  the  design  of  a  lowpass  filter  by  the  windowing 
method.  Use  APL  to  do  this  problem.  Use  STATGRAPHICS  to  generate  3-D 
and  contour  plots. 

1.  Start  with  the  ideal  frequency  response  denned  by  the  following  equation: 

[  0     ,  otherwise  —  re  <  u\ ,u2  <  7r 

Construct  an  APL  array  with  its  values.  Note  that  for  compatibility  with  the 
FFT,  the  range  on  u>i  and  uj2  should  be  zero  to  2  rather  than  — 7r  to  +tt.  Use 
a  32  by  32  point  array. 

2.  Find  the  ideal  impulse  response  2'(nl,n2)  and  plot  it. 

3.  Choose  one  of  the  1-D  windows  and  plot  it.  Define  a  2-D  rectangular  window 
from  it.  Multiply  i(nl,n2)  by  the  window  to  get  /i(nl,n2).  The  window 
should  by  11  by  11  points  in  size. 

4.  Use  the  Fourier  Transform  (FFT)  of  /i(nl,n2)  to  get  the  filter  frequency  re- 
sponse. Plot  this  using  both  3-D  plotting  and  contour  plotting. 


SOLUTION 

1.  By  using  the  function  ILPFILT  (Ideal  LP  filter)  samples  of  the  frequency 
response  of  an  ideal  lowpass  filter  is  generated  (32  by  32  point  array).  The 
plots  achieved  by  using  the  surface  plotting  procedure  in  STATGRAPHICS 
can  be  seen  in  Figure  D.44. 

2.  By  using  the  IFFT2D  function  (Inverse  2-D  Fast  Fourier  Transform)  the  ideal 
impulse  response  is  founded  and  plotted  in  Figure  D.45. 

3.  By  using  the  function  HAMMINGR  we  generate  a  2-D  HAMMING  window 
with  rectangular  base.  See  Figure  D.46. 

4.  Figure  D.47  shows  the  product  of  the  HAMMING  window  and  the  ideal  im- 
pulse response. 
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By  using  the  function  FFT2D  (2-D  Fast  Fourier  Transform)  the  frequency 
response  of  the  windowed  filter  is  found,  and  by  using  the  function  XMAGN 
the  magnitude  of  the  result  is  calculated  and  plotted  as  in  Figures  D.48  and 
D.49  (contour  plot).  Note  the  use  of  the  function  SHIFT2D  to  show  the 
frequency  response  in  the  interval  0  to  2n. 
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Figure  D.44.    Calculating  and  Plotting  the  Ideal  Frequency  Respo 
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Figure  D.45.    Calculating  and  Plotting  the  Ideal  Impulse  Response. 
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Figure  D.46.    Generating  and  Plotting  the  Hamming  Window. 
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Figure  D.47.    The  Product  of  the  Hamming  Window  and  the  Ideal  Im- 
pulse Response. 
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Figure  D.48.    The  Design  Frequency  Response. 
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Figure  D.49.    The  Shifted  Design  Frequency  Response. 
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